From 22263e2adac3a8b95a26ed5cda1ccf353ff558be Mon Sep 17 00:00:00 2001 From: LudwigBoess Date: Wed, 8 Apr 2026 04:14:21 +0000 Subject: [PATCH 1/4] polar coordinates plot --- Project.toml | 8 ++ src/EntityTools.jl | 8 ++ src/plotting/polar.jl | 229 ++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 245 insertions(+) create mode 100644 src/plotting/polar.jl diff --git a/Project.toml b/Project.toml index aaeb2df..d51214a 100644 --- a/Project.toml +++ b/Project.toml @@ -9,6 +9,10 @@ CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0" DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab" Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" FFTW = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341" +GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326" +Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" +LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" +Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7" ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca" Roots = "f2b01f46-fcfa-551c-844a-d8ac1e96c665" @@ -23,6 +27,10 @@ CairoMakie = ">=0.13" DelimitedFiles = "1" Distributions = ">=0.25" FFTW = "1.9.0" +GeometryBasics = "0.5.10" +Interpolations = "0.16.2" +LinearAlgebra = "1.12.0" +Makie = "0.24.9" Printf = "1" ProgressMeter = "1" Roots = "2" diff --git a/src/EntityTools.jl b/src/EntityTools.jl index c4dc7d4..13be0c1 100644 --- a/src/EntityTools.jl +++ b/src/EntityTools.jl @@ -12,6 +12,10 @@ module EntityTools using TOML using FFTW using WriteVTK + using Interpolations + using LinearAlgebra + using Makie + using GeometryBasics include("IO/adios/adios.jl") include("IO/adios/fields.jl") @@ -29,6 +33,7 @@ module EntityTools include("units/time.jl") include("units/Bfield.jl") include("units/temperature.jl") + include("plotting/polar.jl") export phase_map, spectrum, @@ -42,5 +47,8 @@ module EntityTools bp_to_vtk, power_spectrum, get_Temp, + polar_fieldlines!, + polar_pcolor!, + polar_contour! end diff --git a/src/plotting/polar.jl b/src/plotting/polar.jl new file mode 100644 index 0000000..3859192 --- /dev/null +++ b/src/plotting/polar.jl @@ -0,0 +1,229 @@ +""" + dipole_sampling(; nth::Int=30, pole::Float64=1/16) + +Returns an array of angles sampled from a dipole distribution. +""" +function dipole_sampling(; nth::Int=30, pole::Float64=1/16) + nth_poles = floor(Int, nth * pole) + nth_equator = div(nth - 2 * nth_poles, 2) + + # Generate ranges mimicking np.linspace behavior (omitting endpoints where appropriate) + p1 = range(0, stop=π * pole, length=nth_poles + 1)[2:end] + p2 = range(π * pole, stop=π / 2, length=nth_equator + 2)[2:end-1] + p3 = range(π * (1 - pole), stop=π, length=nth_poles + 1)[1:end-1] + + return vcat(p1, p2, p3) +end + +""" + monopole_sampling(; nth::Int=30) + +Returns an array of angles sampled from a monopole distribution. +""" +function monopole_sampling(; nth::Int=30) + return range(0, stop=π, length=nth + 2)[2:end-1] +end + +""" + compute_fieldlines(r_coords, th_coords, fr_data, fth_data, start_points; kwargs...) + +Computes field lines of a vector field defined by radial (`fr`) and azimuthal (`fth`) components. +""" +function compute_fieldlines(r_coords::AbstractVector, th_coords::AbstractVector, + fr_data::AbstractMatrix, fth_data::AbstractMatrix, + start_points::AbstractVector; + direction::String="both", + stop_when::Function=(xy, rth)->false, + ds::Float64=0.1, + maxsteps::Int=1000) + + # Calculate vector field components in Cartesian space for integration + # Note: Assuming data matrices are ordered (th, r) to match original code structure + fxs = zeros(length(th_coords), length(r_coords)) + fys = zeros(length(th_coords), length(r_coords)) + + for (j, r) in enumerate(r_coords) + for (i, th) in enumerate(th_coords) + fxs[i, j] = fr_data[i, j] * sin(th) + fth_data[i, j] * cos(th) + fys[i, j] = fr_data[i, j] * cos(th) - fth_data[i, j] * sin(th) + end + end + + # Nearest neighbor interpolation with 0 fill value + interp_fx = extrapolate(interpolate((th_coords, r_coords), fxs, Gridded(Constant())), 0.0) + interp_fy = extrapolate(interpolate((th_coords, r_coords), fys, Gridded(Constant())), 0.0) + + rmin, rmax = minimum(r_coords), maximum(r_coords) + + # Internal integration step + function _integrate(r_th_start, delta) + r0, th0 = r_th_start + xy = [r0 * sin(th0), r0 * cos(th0)] + rth = [r0, th0] + fieldline = [xy] + + for _ in 1:maxsteps + x, y = xy[1], xy[2] + r = sqrt(x^2 + y^2) + th = atan(-y, x) + π / 2 + rth = [r, th] + + vx = interp_fx(th, r) + vy = interp_fy(th, r) + vmag = sqrt(vx^2 + vy^2) + + if vmag == 0 || isnan(vmag) + break + end + + xy = xy .+ delta .* [vx, vy] ./ vmag + + stop_condition = stop_when(xy, rth) || (rth[1] < rmin) || (rth[1] > rmax) || + (rth[2] < 0) || (rth[2] > π) + + if stop_condition || any(isnan.(xy)) || any(isinf.(xy)) + break + else + push!(fieldline, xy) + end + end + return fieldline + end + + lines = [] + for pt in start_points + if direction == "forward" + push!(lines, _integrate(pt, ds)) + elseif direction == "backward" + push!(lines, _integrate(pt, -ds)) + else + f1 = _integrate(pt, ds) + f2 = _integrate(pt, -ds) + # Combine lines, avoiding duplication of the starting point + push!(lines, vcat(reverse(f2)[1:end-1], f1)) + end + end + + return lines +end + + +""" + polar_fieldlines!(ax, r_coords, th_coords, fr_data, fth_data; kwargs...) + +Plots field lines on an existing Makie axis. +""" +function polar_fieldlines!(ax::Axis, r_coords, th_coords, fr_data, fth_data; + start_points=nothing, + sample_template=nothing, + invert_x::Bool=false, + invert_y::Bool=false, + direction::String="both", + ds::Float64=0.1, + maxsteps::Int=1000, + line_kwargs...) + + if isnothing(start_points) && isnothing(sample_template) + throw(ArgumentError("Either start_points or sample_template must be specified")) + elseif isnothing(start_points) + radius = get(sample_template, :radius, 1.5) + template = get(sample_template, :template, "dipole") + + if template == "dipole" + start_points = [[radius, th] for th in dipole_sampling(; sample_template...)] + elseif template == "monopole" + start_points = [[radius, th] for th in monopole_sampling(; sample_template...)] + else + throw(ArgumentError("Unknown sampling template: $template")) + end + end + + lines_data = compute_fieldlines(r_coords, th_coords, fr_data, fth_data, start_points; + direction=direction, ds=ds, maxsteps=maxsteps) + + for fl in lines_data + fl_matrix = reduce(hcat, fl)' # Convert vector of vectors to matrix + + x_vals = invert_x ? -fl_matrix[:, 1] : fl_matrix[:, 1] + y_vals = invert_y ? -fl_matrix[:, 2] : fl_matrix[:, 2] + + lines!(ax, x_vals, y_vals; line_kwargs...) + end +end + +""" + polar_pcolor!(ax, r_coords, th_coords, data_values; kwargs...) + +Plots a pseudocolor plot of 2D polar data onto a rectilinear projection. +Creates a triangular mesh to cleanly handle non-uniform cells. +""" +function polar_pcolor!(ax::Axis, r_coords::AbstractVector, th_coords::AbstractVector, + data_values::AbstractMatrix; + invert_x::Bool=false, + invert_y::Bool=false, + colormap=:viridis, + mesh_kwargs...) + + nr = length(r_coords) + nth = length(th_coords) + + vertices = Point2f[] + faces = GLTriangleFace[] + colors = Float64[] + + # Map polar to rectilinear coordinates + for (i, th) in enumerate(th_coords) + for (j, r) in enumerate(r_coords) + x = r * sin(th) + y = r * cos(th) + x = invert_x ? -x : x + y = invert_y ? -y : y + push!(vertices, Point2f(x, y)) + push!(colors, data_values[j, i]) + end + end + + # Generate triangulation for the grid + for j in 1:(nth-1) + for i in 1:(nr-1) + idx1 = (j-1)*nr + i + idx2 = (j-1)*nr + i + 1 + idx3 = j*nr + i + 1 + idx4 = j*nr + i + + push!(faces, GLTriangleFace(idx1, idx2, idx3)) + push!(faces, GLTriangleFace(idx1, idx3, idx4)) + end + end + + ax.aspect = DataAspect() + # Use standard Makie mesh rendering for tripped-pseudocolor mapping + mesh!(ax, vertices, faces, color=colors, colormap=colormap; shading=NoShading, mesh_kwargs...) +end + +""" + polar_contour!(ax, r_coords, th_coords, data_values; kwargs...) + +Plots standard contours for 2D polar data on rectilinear axes. +""" +function polar_contour!(ax::Axis, r_coords::AbstractVector, th_coords::AbstractVector, data_values::AbstractMatrix; + invert_x::Bool=false, + invert_y::Bool=false, + contour_kwargs...) + + # Create mapped X and Y grids + X = zeros(length(th_coords), length(r_coords)) + Y = zeros(length(th_coords), length(r_coords)) + + for (j, r) in enumerate(r_coords) + for (i, th) in enumerate(th_coords) + x = r * sin(th) + y = r * cos(th) + X[i, j] = invert_x ? -x : x + Y[i, j] = invert_y ? -y : y + end + end + + ax.aspect = DataAspect() + contour!(ax, X, Y, data_values; contour_kwargs...) +end \ No newline at end of file From 4774e3c05c505087b047e2bc208f45f7f9a85d07 Mon Sep 17 00:00:00 2001 From: LudwigBoess Date: Sat, 16 May 2026 13:56:52 +0000 Subject: [PATCH 2/4] stdout and adios2 performance parsing --- Project.toml | 4 + src/EntityTools.jl | 21 ++- src/IO/adios/performance.jl | 78 ++++++++++ src/IO/parse_output.jl | 283 ++++++++++++++++++++++++++++++++++++ 4 files changed, 385 insertions(+), 1 deletion(-) create mode 100644 src/IO/adios/performance.jl create mode 100644 src/IO/parse_output.jl diff --git a/Project.toml b/Project.toml index d51214a..5a471ec 100644 --- a/Project.toml +++ b/Project.toml @@ -11,11 +11,13 @@ Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" FFTW = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341" GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326" Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" +JSON3 = "0f8b85d8-7281-11e9-16c2-39a750bddbf1" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7" ProgressMeter = "92933f4c-e287-5a05-a399-4b506db050ca" Roots = "f2b01f46-fcfa-551c-844a-d8ac1e96c665" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" TOML = "fa267f1f-6049-4f14-aa54-33bafae1ed76" Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" @@ -29,11 +31,13 @@ Distributions = ">=0.25" FFTW = "1.9.0" GeometryBasics = "0.5.10" Interpolations = "0.16.2" +JSON3 = "1.14.3" LinearAlgebra = "1.12.0" Makie = "0.24.9" Printf = "1" ProgressMeter = "1" Roots = "2" +Statistics = "1" StatsBase = ">=0.34" TOML = "1.0.3" Unitful = "1" diff --git a/src/EntityTools.jl b/src/EntityTools.jl index 13be0c1..89eb6d3 100644 --- a/src/EntityTools.jl +++ b/src/EntityTools.jl @@ -3,6 +3,7 @@ module EntityTools using CairoMakie using Printf using StatsBase + using Statistics using Unitful using Distributions using StatsBase @@ -11,17 +12,20 @@ module EntityTools using Printf using TOML using FFTW - using WriteVTK + using WriteVTK using Interpolations using LinearAlgebra using Makie using GeometryBasics + using JSON3 include("IO/adios/adios.jl") include("IO/adios/fields.jl") include("IO/adios/particles.jl") include("IO/adios/spectra.jl") + include("IO/adios/performance.jl") include("IO/vtk/write_vtk.jl") + include("IO/parse_output.jl") include("calc/phase.jl") include("calc/spectra.jl") include("calc/powerspectrum.jl") @@ -38,6 +42,21 @@ module EntityTools export phase_map, spectrum, parse_timing, + parse_output, + SimOutput, + StepData, + substep_series, + total_series, + active_particles, + timestep_durations, + elapsed_times, + remaining_times, + sim_times, + step_numbers, + species_series, + summary_table, + parse_compound_time, + adios2_throughput, EntityData, EntityUnits, find_closest_time, diff --git a/src/IO/adios/performance.jl b/src/IO/adios/performance.jl new file mode 100644 index 0000000..491e926 --- /dev/null +++ b/src/IO/adios/performance.jl @@ -0,0 +1,78 @@ +""" + adios2_throughput(path::AbstractString) + +Read an ADIOS2 BP5 `profiling.json` and return a NamedTuple with: + +- `aggregate_TBps` : Σwbytes(all ranks) / max(durable_us) [decimal TB/s] +- `per_rank_mean_GBps` : mean over writer ranks of wbytes / durable_us [decimal GB/s] +- `per_rank_std_GBps` : stddev of the same +- `per_rank_raw_GBps` : mean of wbytes / Σwrite.mus over writer ranks [decimal GB/s] — + the raw POSIX write() rate, kept for reference; this measures + the rate of memcpy into the page/client cache and overstates + durable bandwidth, often above the filesystem spec +- `n_writers` : number of ranks that issued writes (BP5 aggregators) +- `n_ranks` : total ranks in the profile + +`durable_us` per rank is `ES + ES_close + DC_WaitOnAsync1 + DC_WaitOnAsync2`, +i.e. the EndStep window plus the Close-time drain of BP5's async writer thread. +With `AsyncWrite=ON`, `EndStep` returns before bytes are durable; the drain shows +up as `DC_WaitOnAsync*` at Close. Using this window gives bandwidth bounded by +what actually reached storage. + +Aggregate throughput uses `max(durable_us)` over all ranks (writers and +non-writers) — non-writer ranks still block at Close waiting for the aggregator, +and the wall clock is bounded by the slowest rank. +""" +function adios2_throughput(path::AbstractString) + raw = read(path, String) + # ADIOS2 sometimes emits a trailing comma before the closing array bracket + raw = replace(raw, r",(\s*\])" => s"\1") + data = JSON3.read(raw) + + total_bytes = 0 + rank_eff_GBps = Float64[] + rank_raw_GBps = Float64[] + durable_times_us = Float64[] + + for rec in data + rank_bytes = 0 + rank_write_us = 0 + for (k, v) in pairs(rec) + startswith(String(k), "transport_") || continue + v isa JSON3.Object || continue + rank_bytes += get(v, :wbytes, 0) + wr = get(v, :write, nothing) + wr === nothing || (rank_write_us += get(wr, :mus, 0)) + end + total_bytes += rank_bytes + + es = Float64(get(rec, :ES_mus, 0)) + esclose = Float64(get(rec, :ES_close_mus, 0)) + wait1 = Float64(get(rec, :DC_WaitOnAsync1_mus, 0)) + wait2 = Float64(get(rec, :DC_WaitOnAsync2_mus, 0)) + durable_us = es + esclose + wait1 + wait2 + durable_us > 0 && push!(durable_times_us, durable_us) + + if rank_bytes > 0 && durable_us > 0 + # bytes / μs == MB/s (decimal); /1000 -> GB/s + push!(rank_eff_GBps, rank_bytes / durable_us / 1000) + end + if rank_bytes > 0 && rank_write_us > 0 + push!(rank_raw_GBps, rank_bytes / rank_write_us / 1000) + end + end + + isempty(durable_times_us) && error("no ES/Close/Async timing entries found in $path") + wall_us = maximum(durable_times_us) + # bytes / μs == MB/s; /1e6 -> TB/s (decimal) + aggregate_TBps = total_bytes / wall_us / 1e6 + + return ( + aggregate_TBps = aggregate_TBps, + per_rank_mean_GBps = isempty(rank_eff_GBps) ? NaN : mean(rank_eff_GBps), + per_rank_std_GBps = length(rank_eff_GBps) < 2 ? NaN : std(rank_eff_GBps), + per_rank_raw_GBps = isempty(rank_raw_GBps) ? NaN : mean(rank_raw_GBps), + n_writers = length(rank_eff_GBps), + n_ranks = length(data), + ) +end diff --git a/src/IO/parse_output.jl b/src/IO/parse_output.jl new file mode 100644 index 0000000..a948ca1 --- /dev/null +++ b/src/IO/parse_output.jl @@ -0,0 +1,283 @@ +# --------------------------------------------------------------------------- +# Time-unit handling +# --------------------------------------------------------------------------- + +const TIME_UNITS = Dict( + "ns" => 1e-9, + "µs" => 1e-6, + "us" => 1e-6, + "ms" => 1e-3, + "s" => 1.0, + "min" => 60.0, + "h" => 3600.0, +) + +to_seconds(value::Real, unit::AbstractString) = value * TIME_UNITS[unit] + +const TIME_PAIR = r"([\d.]+)\s*(ns|µs|us|ms|min|h|s)" + +""" + parse_compound_time(str) + +Parse compound durations like `"10s 143ms"`, `"3min 39s"`, `"286ms 625µs"` +into seconds. +""" +function parse_compound_time(str::AbstractString) + total = 0.0 + for m in eachmatch(TIME_PAIR, str) + total += to_seconds(parse(Float64, m.captures[1]), m.captures[2]) + end + return total +end + +# --------------------------------------------------------------------------- +# Data containers +# --------------------------------------------------------------------------- + +""" + StepData + +Per-step record from an entity `.out` file. All times are in seconds. +Substep timings are accessible by name: `step["CurrentDeposit"]`. +""" +struct StepData + step::Int + sim_time::Float64 + dt::Float64 + species::Dict{String,Float64} # label => particle count (global total) + n_active::Float64 # Σ species counts + substeps::Dict{String,Float64} # substep name => seconds + pre_total_substeps::Vector{String} # order of substeps included in Total + post_total_substeps::Vector{String} # order of substeps reported after Total + total::Float64 # "Total" line in seconds + timestep_duration::Float64 # seconds + remaining_time::Float64 # seconds + elapsed_time::Float64 # seconds +end + +Base.getindex(s::StepData, name::AbstractString) = s.substeps[name] +Base.haskey(s::StepData, name::AbstractString) = haskey(s.substeps, name) + +""" + SimOutput + +Container for all per-step records parsed from a `.out` file. Iterable and +indexable by step index. The set of substep and species names observed across +the whole run is kept in `substep_names` / `species_labels`. +""" +struct SimOutput + steps::Vector{StepData} + substep_names::Vector{String} # union of substep names, insertion order + species_labels::Vector{String} # union of species labels, insertion order + source::String # path the data was parsed from +end + +Base.length(s::SimOutput) = length(s.steps) +Base.getindex(s::SimOutput, i::Integer) = s.steps[i] +Base.iterate(s::SimOutput, st...) = iterate(s.steps, st...) +Base.firstindex(s::SimOutput) = 1 +Base.lastindex(s::SimOutput) = length(s.steps) + +# --------------------------------------------------------------------------- +# Parsing +# --------------------------------------------------------------------------- + +const STEP_LINE = r"^Step:\s*(\d+)" +const NUMBER = raw"[0-9]+(?:\.[0-9]+)?(?:[eE][+-]?[0-9]+)?" +const TIME_LINE = Regex("^Time:\\s*($NUMBER).*?Δt\\s*=\\s*($NUMBER)") +const SUBSTEP_LN = r"^\s+([A-Za-z][A-Za-z0-9_]*)\.{2,}\s*([\d.]+)\s*(ns|µs|us|ms|min|h|s)\b" +const TOTAL_LINE = r"^Total\s+([\d.]+)\s*(ns|µs|us|ms|min|h|s)\b" +const SPECIES_LN = r"^\s+species\s+\d+\s+\(([^)]+)\)\.{2,}\s*([\d.eE+-]+)" +const TIMESTEPDUR = r"^Timestep duration:\s*(.+)$" +const REMAINING = r"^Remaining time:\s*(.+)$" +const ELAPSED = r"^Elapsed time:\s*(.+)$" + +""" + parse_output(path) -> SimOutput + +Parse an entity `.out` log file. The number and names of substeps are +discovered dynamically. All times are stored in seconds. +""" +function parse_output(path::AbstractString) + isfile(path) || error("file not found: $path") + + steps = StepData[] + substep_order = String[] + species_order = String[] + + in_step = false + cur_step = 0 + cur_simtime = 0.0 + cur_dt = 0.0 + cur_species = Dict{String,Float64}() + cur_substeps = Dict{String,Float64}() + cur_pre = String[] + cur_post = String[] + cur_total = NaN + cur_tsdur = NaN + cur_remaining = NaN + cur_elapsed = NaN + saw_total = false + + flush_step! = function () + n_active = isempty(cur_species) ? 0.0 : sum(values(cur_species)) + push!(steps, StepData( + cur_step, cur_simtime, cur_dt, + copy(cur_species), n_active, + copy(cur_substeps), copy(cur_pre), copy(cur_post), + cur_total, cur_tsdur, cur_remaining, cur_elapsed, + )) + end + + reset_step! = function () + empty!(cur_species); empty!(cur_substeps) + empty!(cur_pre); empty!(cur_post) + cur_total = NaN; cur_tsdur = NaN + cur_remaining = NaN; cur_elapsed = NaN + saw_total = false + end + + open(path, "r") do io + for line in eachline(io) + m = match(STEP_LINE, line) + if m !== nothing + in_step && flush_step!() + reset_step!() + in_step = true + cur_step = parse(Int, m.captures[1]) + continue + end + in_step || continue + + m = match(TIME_LINE, line) + if m !== nothing + cur_simtime = parse(Float64, m.captures[1]) + cur_dt = parse(Float64, m.captures[2]) + continue + end + + m = match(TOTAL_LINE, line) + if m !== nothing + cur_total = to_seconds(parse(Float64, m.captures[1]), m.captures[2]) + saw_total = true + continue + end + + m = match(SUBSTEP_LN, line) + if m !== nothing + name = m.captures[1] + t = to_seconds(parse(Float64, m.captures[2]), m.captures[3]) + cur_substeps[name] = t + if saw_total + push!(cur_post, name) + else + push!(cur_pre, name) + end + name in substep_order || push!(substep_order, name) + continue + end + + m = match(SPECIES_LN, line) + if m !== nothing + label = m.captures[1] + cur_species[label] = parse(Float64, m.captures[2]) + label in species_order || push!(species_order, label) + continue + end + + m = match(TIMESTEPDUR, line) + if m !== nothing + cur_tsdur = parse_compound_time(m.captures[1]) + continue + end + + m = match(REMAINING, line) + if m !== nothing + cur_remaining = parse_compound_time(m.captures[1]) + continue + end + + m = match(ELAPSED, line) + if m !== nothing + cur_elapsed = parse_compound_time(m.captures[1]) + continue + end + end + in_step && flush_step!() + end + + return SimOutput(steps, substep_order, species_order, abspath(path)) +end + +# --------------------------------------------------------------------------- +# Series accessors +# --------------------------------------------------------------------------- + +""" + substep_series(sim, name) + +Per-step time series (seconds) for substep `name`. Steps where the substep +is absent return `NaN`. +""" +substep_series(sim::SimOutput, name::AbstractString) = + [get(s.substeps, name, NaN) for s in sim.steps] + +total_series(sim::SimOutput) = [s.total for s in sim.steps] +active_particles(sim::SimOutput) = [s.n_active for s in sim.steps] +timestep_durations(sim::SimOutput) = [s.timestep_duration for s in sim.steps] +elapsed_times(sim::SimOutput) = [s.elapsed_time for s in sim.steps] +remaining_times(sim::SimOutput) = [s.remaining_time for s in sim.steps] +sim_times(sim::SimOutput) = [s.sim_time for s in sim.steps] +step_numbers(sim::SimOutput) = [s.step for s in sim.steps] + +species_series(sim::SimOutput, label::AbstractString) = + [get(s.species, label, NaN) for s in sim.steps] + +# --------------------------------------------------------------------------- +# Aggregates +# --------------------------------------------------------------------------- + +_clean(xs) = filter(!isnan, xs) +_mean(xs) = isempty(xs) ? NaN : sum(xs) / length(xs) + +""" + sum(sim::SimOutput) -> seconds + sum(sim::SimOutput, name) -> seconds + maximum(sim::SimOutput, name) -> seconds + minimum(sim::SimOutput, name) -> seconds + mean(sim::SimOutput) -> seconds + mean(sim::SimOutput, name) -> seconds + +Aggregate over the simulation. With no `name`, the aggregate runs over the +per-step `Total` line. With a `name`, it runs over the named substep +(e.g. `"CurrentDeposit"`); steps where the substep is absent are skipped. +""" +Base.sum(sim::SimOutput) = sum(_clean(total_series(sim)); init = 0.0) + +Base.sum(sim::SimOutput, name::AbstractString) = + sum(_clean(substep_series(sim, name)); init = 0.0) + +Base.maximum(sim::SimOutput, name::AbstractString) = + maximum(_clean(substep_series(sim, name)); init = -Inf) + +Base.minimum(sim::SimOutput, name::AbstractString) = + minimum(_clean(substep_series(sim, name)); init = Inf) + +Statistics.mean(sim::SimOutput) = _mean(_clean(total_series(sim))) + +Statistics.mean(sim::SimOutput, name::AbstractString) = + _mean(_clean(substep_series(sim, name))) + +""" + summary_table(sim) -> Vector{NamedTuple} + +One row per substep with `(sum, mean, min, max)` in seconds, in the order the +substeps first appeared in the file. +""" +function summary_table(sim::SimOutput) + [(name = n, + sum = sum(sim, n), + mean = Statistics.mean(sim, n), + min = minimum(sim, n), + max = maximum(sim, n)) for n in sim.substep_names] +end From 66c18b6d70269e35a6c7273dce49ddafbb022ed9 Mon Sep 17 00:00:00 2001 From: LudwigBoess Date: Sat, 16 May 2026 13:57:01 +0000 Subject: [PATCH 3/4] some initial tests --- test/runtests.jl | 74 ++++++++++++++++++++++++++++++++++++++++++++++-- 1 file changed, 71 insertions(+), 3 deletions(-) diff --git a/test/runtests.jl b/test/runtests.jl index 3f76732..4afd840 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,6 +1,74 @@ -using EntityTools using Test +using EntityTools # Replace with your package's actual module name if different +using Makie @testset "EntityTools.jl" begin - # Write your tests here. -end + + @testset "Sampling Methods" begin + @testset "Monopole Sampling" begin + nth = 30 + mono = monopole_sampling(nth=nth) + + @test length(mono) == nth + @test all(mono .> 0.0) + @test all(mono .< π) + @test issorted(mono) + end + + @testset "Dipole Sampling" begin + nth = 30 + pole = 1/16 + dip = dipole_sampling(nth=nth, pole=pole) + + # Math logic from the Python translation + expected_poles = floor(Int, nth * pole) + expected_equator = div(nth - 2 * expected_poles, 2) + expected_length = 2 * expected_poles + expected_equator + + @test length(dip) == expected_length + @test all(dip .> 0.0) + @test all(dip .< π) + @test issorted(dip) + end + end + + @testset "Fieldline Integration" begin + # Create a simple mock grid and uniform vector field + r_coords = collect(1.0:0.5:5.0) + th_coords = collect(0.0:0.2:π) + + fr_data = ones(length(th_coords), length(r_coords)) + fth_data = zeros(length(th_coords), length(r_coords)) + + start_pts = [[2.0, π/4], [3.0, π/2]] + + lines = compute_fieldlines(r_coords, th_coords, fr_data, fth_data, start_pts, maxsteps=10) + + @test length(lines) == 2 + @test length(lines[1]) > 1 + @test length(lines[1][1]) == 2 # Output should be [x, y] arrays + @test typeof(lines) <: AbstractVector + end + + @testset "Plotting Recipes (Makie)" begin + # Test that plotting functions execute without errors + fig = Figure() + ax = Axis(fig[1, 1]) + + r_coords = collect(1.0:1.0:5.0) + th_coords = collect(0.0:0.5:π) + data = rand(length(th_coords), length(r_coords)) + + fr_data = rand(length(th_coords), length(r_coords)) + fth_data = rand(length(th_coords), length(r_coords)) + + @test_nowarn polar_pcolor!(ax, r_coords, th_coords, data) + @test_nowarn polar_contour!(ax, r_coords, th_coords, data) + + # Test fieldlines plot with the dictionary template + sample_dict = Dict(:template => "monopole", :nth => 5, :radius => 2.0) + @test_nowarn plot_fieldlines!(ax, r_coords, th_coords, fr_data, fth_data; + sample_template=sample_dict, maxsteps=5) + end + +end \ No newline at end of file From 6ccfb09dbf1b69014efc88357456d17155849173 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ludwig=20B=C3=B6ss?= Date: Sun, 17 May 2026 16:54:06 -0500 Subject: [PATCH 4/4] softer LinearAlgebra requirements --- Project.toml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index 5a471ec..6107473 100644 --- a/Project.toml +++ b/Project.toml @@ -32,7 +32,7 @@ FFTW = "1.9.0" GeometryBasics = "0.5.10" Interpolations = "0.16.2" JSON3 = "1.14.3" -LinearAlgebra = "1.12.0" +LinearAlgebra = "1" Makie = "0.24.9" Printf = "1" ProgressMeter = "1" @@ -41,7 +41,7 @@ Statistics = "1" StatsBase = ">=0.34" TOML = "1.0.3" Unitful = "1" -WriteVTK = "1.21.2" +WriteVTK = "1" julia = ">=1.10" [extras]