Tarang.jl

Features

Quick Start

Installation

using Pkg
Pkg.add(url="https://github.com/subhk/Tarang.jl")

That is the whole installation. MPI, PencilArrays, PencilFFTs and KernelAbstractions are hard dependencies — they are installed with Tarang, and distributed runs work out of the box.

GPU support is the one opt-in: CUDA is a weak dependency loaded through a package extension, so add it only if you want it.

Pkg.add("CUDA")     # enables the GPU backend (TarangCUDAExt)

1D Diffusion

using Tarang

domain = PeriodicDomain(64)                     # 64-point periodic [0, 2π]
T = ScalarField(domain, "T")                    # Temperature field

problem = IVP([T])
add_parameters!(problem, kappa=0.01)
add_equation!(problem, "∂t(T) - kappa*Δ(T) = 0")

set!(T, x -> sin(x))                            # Initial condition
solver = InitialValueSolver(problem, RK222(); dt=0.01)
run!(solver; stop_time=1.0)                     # That's it!

The single mode sin(x) decays as exp(-κt); after t = 1 the solver gives max|T| = 0.99004983 against the exact 0.99004983 — a relative error of 4e-12.

2D Rayleigh-Benard Convection

A bounded (Chebyshev) direction needs the tau method: one tau variable per boundary condition, lifted into the equations. See The Tau Method for Boundary Conditions for the why.

using Tarang

Lx, Lz = 4.0, 1.0
Nx, Nz = 64, 32
Rayleigh, Prandtl = 2e4, 1.0

coords = CartesianCoordinates("x", "z")
dist   = Distributor(coords; dtype=Float64, device=CPU())
xbasis = RealFourier(coords["x"]; size=Nx, bounds=(0.0, Lx), dealias=3/2)
zbasis = ChebyshevT(coords["z"]; size=Nz, bounds=(0.0, Lz), dealias=3/2)
domain = Domain(dist, (xbasis, zbasis))

p = ScalarField(domain, "p")                   # Pressure
T = ScalarField(domain, "T")                   # Temperature
u = VectorField(domain, "u")                   # Velocity

# One tau per boundary condition, carrying the Fourier bases (one per x-mode).
tau_p  = ScalarField(dist, "tau_p",  (), Float64)          # pressure gauge
tau_T1 = ScalarField(dist, "tau_T1", (xbasis,), Float64)
tau_T2 = ScalarField(dist, "tau_T2", (xbasis,), Float64)
tau_u1 = VectorField(dist, coords, "tau_u1", (xbasis,), Float64)
tau_u2 = VectorField(dist, coords, "tau_u2", (xbasis,), Float64)

# First-order reduction: grad_X = ∇X + ẑ·lift(τ)
ex, ez     = unit_vector_fields(coords, dist)
lift_basis = derivative_basis(zbasis, 1)
τ_lift(A)  = lift(A, lift_basis, -1)
grad_u = grad(u) + ez * τ_lift(tau_u1)
grad_T = grad(T) + ez * τ_lift(tau_T1)

problem = IVP([p, T, u, tau_p, tau_T1, tau_T2, tau_u1, tau_u2])
add_parameters!(problem, nu=Prandtl, buoy=Rayleigh*Prandtl, ez=ez,
                grad_u=grad_u, grad_T=grad_T, τ_lift=τ_lift)

add_equation!(problem, "trace(grad_u) + tau_p = 0")                                  # continuity
add_equation!(problem, "∂t(T) - div(grad_T) + τ_lift(tau_T2) = -u⋅∇(T)")
add_equation!(problem, "∂t(u) - nu*div(grad_u) + ∇(p) - buoy*T*ez + τ_lift(tau_u2) = -u⋅∇(u)")

# Boundary conditions go through add_bc!, never add_equation!.
add_bc!(problem, "T(z=0) = 1")                 # hot bottom
add_bc!(problem, "T(z=$Lz) = 0")               # cold top
add_bc!(problem, "u(z=0) = 0")                 # no-slip
add_bc!(problem, "u(z=$Lz) = 0")
add_bc!(problem, "integ(p) = 0")               # pressure gauge

# Conduction profile 1 - z, plus wall-damped noise to seed convection.
x, z = local_grids(dist, xbasis, zbasis)
fill_random!(T, "g"; seed=42, distribution="normal", scale=1e-3)
get_grid_data(T) .*= z' .* (1.0 .- z')
get_grid_data(T) .+= 1.0 .- z'
ensure_layout!(T, :c)

solver = InitialValueSolver(problem, RK222(); dt=1e-4)
diagnose(solver)                               # Print solver summary
run!(solver; stop_time=0.1, log_interval=100)

Time is measured in thermal diffusion times, so buoy = Ra·Pr and nu = Pr. Over those first 1000 steps the seeded noise grows into rolls (max|u_z| reaches 3.76) while the walls stay pinned to machine precision: max|T(z=0) − 1| = 3.6e-15, max|T(z=Lz)| = 8.1e-15.

GPU Example

using Tarang, CUDA

# Just add device=GPU() — everything else stays the same
domain = PeriodicDomain(512, 512; device=GPU(), dtype=Float32)
field = ScalarField(domain, "u")
forward_transform!(field)   # Uses cuFFT automatically

device is the only line that changes: the same script with device=CPU() runs on the CPU, and dtype=Float32 gives ComplexF32 coefficients either way. GPU() needs CUDA.jl loaded — without it the constructor throws rather than silently falling back.

Scientific Applications

Problem Types

Spectral Bases

Documentation

• Installation
• Tutorials Overview
• GPU Computing
• Parallelism
• Coordinates API
• Example Gallery

Installation Options

MPI, MPIPreferences, PencilArrays, PencilFFTs and KernelAbstractions are [deps] of Tarang, so a plain install is already MPI-capable: mpiexec -n 4 julia --project=. run.jl works with the MPI binary that MPI.jl ships. Only CUDA is a [weakdeps] package extension.

Contributing

We welcome contributions! See our GitHub repository for:

• Bug reports and feature requests
• Documentation improvements
• Pull requests

Citation

If you use Tarang.jl in your research, please cite:

@software{tarang_jl,
  author = {Kar, Subhajit},
  title  = {Tarang.jl: A Spectral PDE Solver for Julia},
  url    = {https://github.com/subhk/Tarang.jl},
  year   = {2024}
}

License

Tarang.jl is released under the MIT License.