Julia Scientific Computing skill for Claude Code

Done well means: hot functions are type-stable and allocation-free in their inner loops, proved with @code_warntype and @btime rather than assumed; every project runs from its own Project.toml with a committed Manifest.toml for applications; code is organized around functions and dispatch, not scripts mutating globals; and the first-run latency is known and, where it matters, paid once by precompilation.

Versions

• Install and switch Julia with juliaup: juliaup add release, juliaup default release, juliaup status. Stable in late 2026 is 1.13.x; the LTS channel (juliaup add lts) has been 1.10. Pin the minimum supported version in [compat] for packages.
• Start every session in the project: julia --project=. -t auto (-t auto enables threads; Julia starts single-threaded by default).

Ground rules

• Put performance-critical code in functions. Top-level script code and non-const globals are type-unstable and slow.
• Never benchmark the first call; it includes compilation. Use BenchmarkTools with $ interpolation for arguments.
• Struct fields must have concrete types (or type parameters). x::Real or an untyped field makes every access dynamic.
• Do not add methods to functions you do not own for types you do not own (type piracy). It changes behavior for every package loaded in the session.
• Arrays are column-major: the first index varies fastest. Loop for j in axes(A,2), i in axes(A,1).
• Prove before using @inbounds or @fastmath; both trade safety or IEEE semantics for speed. Run tests with --check-bounds=yes regardless.

Workflow

- [ ] 1. Set up the environment (Project.toml, Manifest.toml, Revise)
- [ ] 2. Write it plainly in functions, with tests
- [ ] 3. Measure: @btime / @benchmark on realistic sizes
- [ ] 4. Check type stability: @code_warntype, JET @report_opt
- [ ] 5. Remove allocations in hot loops; profile the rest
- [ ] 6. Parallelize only after the serial version is efficient
- [ ] 7. Reduce latency (precompile workloads, sysimage) if users feel it

1. Environment

julia> ]                      # Pkg REPL
(@v1.13) pkg> activate .      # creates/uses ./Project.toml
(proj) pkg> add DataFrames CSV CairoMakie BenchmarkTools
(proj) pkg> instantiate       # on a fresh clone: install exactly what Manifest.toml records
(proj) pkg> status --outdated

Project type	Commit Manifest.toml?	Compat
Analysis, paper, application, service	Yes (exact reproducibility)	Optional
Registered package / library	No	Required [compat] for every dependency and julia

• Keep the global environment (@v1.x) to dev tools only: Revise, BenchmarkTools, JET, maybe OhMyREPL. Everything else goes in the project.
• using Revise before loading your code so edits reload without restarting.
• Structure non-trivial work as a package (pkg> generate MyAnalysis) with src/ and test/runtests.jl; scripts then just using MyAnalysis.

3. Measure

using BenchmarkTools
x = rand(10^6)
@btime mysum($x)              # min time and allocations; $ avoids global-variable overhead
@benchmark mysum($x)          # distribution (median, GC time)
@allocated mysum(x)           # bytes for one call (after warm-up)

Targets: an inner kernel over plain arrays should report 0 allocations. A loop that allocates per iteration is the most common 10-100x slowdown.

4. Type stability

@code_warntype f(x)           # red Any / Union{...} in Body or variables = instability
using JET
@report_opt f(x)              # lists runtime dispatch sites through the call tree

Common causes and fixes are in references/performance-patterns.md; read it when @code_warntype shows Any or @btime shows unexpected allocations.

Quick rules:

• A variable must keep one type: s = 0 then s += 0.5 makes s a Union. Use s = zero(eltype(x)).
• Return the same type on every branch (nothing vs a number is a Union; small Unions are fine, Any is not).
• Containers: Float64[], not [] (which is Vector{Any}).
• Closures capturing a variable that is reassigned get boxed (Core.Box); assign once or use let.
• Function barriers: when data arrives with an unknown type (parsed file, Dict{String,Any}), pass it to a separate kernel function; dispatch specializes the kernel once.

5. Allocations and profiling

• Preallocate outputs and write mutating versions: mul!(C, A, B), f!(out, x).
• Fuse broadcasts with dots: @. y = a * x + b makes one loop, no temporaries.
• Slices copy: A[:, j] allocates; use @views or view(A, :, j).
• Small fixed-size vectors (3D points, 4x4 matrices): StaticArrays SVector/SMatrix.
• Profile: using Profile; @profile run(); Profile.print(format=:flat, sortedby=:count), or @profview in the VS Code Julia extension for a flame graph.

6. Parallelism

Work shape	Tool
Independent loop iterations, shared memory	Threads.@threads or OhMyThreads.tmap with -t auto
Task-based or uneven work	Threads.@spawn + fetch
Multiple machines or processes	Distributed (addprocs, pmap)
Dense linear algebra	Already multithreaded via BLAS; check BLAS.get_num_threads() and avoid oversubscription with Julia threads
GPU	CUDA.jl / AMDGPU.jl / Metal.jl with broadcast or KernelAbstractions

Never write to shared accumulators from threads without a reduction or per-task buffers; Threads.threadid()-indexed buffers are unsafe because tasks can migrate.

7. Latency (time to first X)

Problem	Fix
Your package’s first call is slow for users	PrecompileTools.@compile_workload with a representative workload inside the package module
Your own scripts load many packages slowly	A local Startup package that usings them and runs a small @compile_workload
CLI tool or deployment must start in well under a second	PackageCompiler create_sysimage with a precompile execution file, or experimental juliac --trim (1.12+) for small static-dispatch programs
Iterating interactively	Keep one session alive with Revise; do not restart

using PrecompileTools: @setup_workload, @compile_workload
@setup_workload begin
    data = rand(100, 3)
    @compile_workload begin
        fit_model(data)
    end
end

Invalidations (using SnoopCompileCore; invalidations = @snoop_invalidations using X) undo precompilation; type piracy and overly generic methods are the usual cause.

Dispatch design

• Define behavior on abstract types and functions, data in concrete structs: abstract type Kernel end; struct Gaussian <: Kernel; σ::Float64; end; evaluate(k::Gaussian, r) = exp(-r^2 / (2k.σ^2)).
• Parametric structs for generic numerics: struct Particle{T<:AbstractFloat}; x::T; v::T; end.
• Accept abstract argument types (AbstractVector{<:Real}) for generality; type annotations on arguments do not make code faster, concrete field and container types do.
• Extend Base generics (Base.length, Base.getindex, Base.show) for your own types to plug into the ecosystem.

DataFrames and plotting

using CSV, DataFrames, Statistics
df = CSV.read("runs.csv", DataFrame; types=Dict(:id => String), missingstring="NA")
clean = dropmissing(df, :energy)
summary = combine(groupby(clean, :material), :energy => mean => :mean_E, nrow => :n)
transform!(clean, [:mass, :v] => ByRow((m, v) -> 0.5m * v^2) => :kinetic)

• Use select/transform/combine with source => function => target; row loops over a DataFrame are type-unstable because columns have heterogeneous types.
• For a hot loop over columns, extract them first (x = df.x) and pass to a kernel.

Plotting need	Choose
Publication figures, PDF/SVG/PNG	CairoMakie
Interactive 3D or large data	GLMakie (needs OpenGL)
Quick plots in notebooks with many backends	Plots.jl

Save explicitly: save("fig.pdf", fig) (Makie) or savefig("fig.pdf") (Plots).

Interop

using PythonCall                      # Python from Julia; deps recorded in CondaPkg.toml
np = pyimport("numpy")
v = pyconvert(Vector{Float64}, np.linspace(0, 1, 5))
a = PyArray(np.zeros(10))             # no-copy view of a NumPy buffer

• Add Python packages with pkg> conda add scipy (after using CondaPkg) and commit CondaPkg.toml. To use an existing interpreter instead, set JULIA_CONDAPKG_BACKEND=Null and JULIA_PYTHONCALL_EXE=/path/to/python.
• From Python, juliacall exposes Julia (from juliacall import Main as jl).
• C and Fortran: @ccall libm.cos(1.0::Cdouble)::Cdouble; pass arrays as Ptr{Cdouble} and keep them rooted with GC.@preserve while C holds the pointer.

Julia vs Python + NumPy

Situation	Better choice
Custom loops that do not vectorize (ODE right-hand sides, Monte Carlo with branching, agent models, custom kernels)	Julia: loops compile to native code
Differential equations, optimization, AD-heavy models	Julia (DifferentialEquations.jl, JuMP, Enzyme/ForwardDiff) is a strong default
Workload is already vectorized NumPy/BLAS calls	Roughly equal; switching gains little
Deep learning training, cloud ML services, broad library needs	Python (PyTorch/JAX ecosystem)
Short scripts run once, cold start matters, team only knows Python	Python; Julia’s compile latency dominates
Mixed team	Python front end calling Julia kernels via juliacall, or the reverse

Do not claim speedups without a benchmark of both versions on the user’s data size.

Gotchas

• Integer overflow wraps silently: 2^64 == 0. Use big(2)^64 or Int128 when ranges can exceed Int64.
• 1:n with n = 0 is empty, and arrays are 1-based; eachindex(A) and axes(A, d) are correct for offset arrays too.
• a = b for arrays aliases; use copy(b).
• == vs isequal/===: NaN == NaN is false; missing == 1 is missing.
• Precompilation can fail silently in a broken environment; run pkg> precompile and read the errors.
• Threads.@threads with a shared push! is a data race; it may “work” on small inputs.

Verification

1. pkg> test passes from a fresh clone after instantiate.
2. @code_warntype on each hot function shows no Any; @report_opt reports nothing in your code.
3. @btime before and after, at the target data size, with allocations; report both.
4. Numerical results match the previous implementation within a stated tolerance (isapprox(a, b; rtol=1e-10)), including on edge-case inputs.

If speed did not improve, profile before guessing; the bottleneck is often I/O, a DataFrame row loop, or a type-unstable function barrier one level up.