Architecture
This document describes Prismo’s architecture, design patterns, and code organization.
Overview
Prismo follows a modular architecture with clear separation of concerns:
┌────────────────────────────────────────────────┐
│ User Interface │
│ (Simulation, convenience APIs) │
└────────────────────────────────────────────────┘
│
┌────────────────────────────────────────────────┐
│ Core FDTD Engine │
│ (Solver, Grid, Fields, Time-stepping) │
└────────────────────────────────────────────────┘
│
┌─────────────┬───────────────┬──────────────────┐
│ Sources │ Monitors │ Boundaries │
│ (Excite) │ (Measure) │ (Absorb/BC) │
└─────────────┴───────────────┴──────────────────┘
│
┌─────────────┬───────────────┬──────────────────┐
│ Materials │ Geometry │ Backends │
│ (ε, μ, σ) │ (Shapes) │ (CPU/GPU) │
└─────────────┴───────────────┴──────────────────┘
Directory Structure
prismo/
├── core/ # Core FDTD engine
│ ├── fields.py # Field storage and access
│ ├── grid.py # Yee grid implementation
│ ├── solver.py # Main FDTD solver
│ └── simulation.py # High-level simulation interface
├── sources/ # Field sources
│ ├── base.py # Abstract source base class
│ ├── point.py # Point sources
│ ├── plane_wave.py # Plane wave sources
│ ├── mode.py # Mode sources
│ └── waveform.py # Temporal waveforms
├── monitors/ # Field monitors
│ ├── base.py # Abstract monitor base class
│ ├── field.py # Field monitors
│ ├── flux.py # Flux/power monitors
│ └── mode_monitor.py # Mode expansion monitors
├── boundaries/ # Boundary conditions
│ ├── pml.py # PML absorbing boundaries
│ └── mode_port.py # Mode port boundaries
├── materials/ # Material models
│ ├── dispersion.py # Dispersive materials
│ ├── library.py # Material library
│ └── tensor.py # Anisotropic materials
├── modes/ # Mode solver
│ └── solver.py # Eigenmode solver
├── geometry/ # Geometric primitives
│ └── shapes.py # Box, sphere, cylinder, etc.
├── backends/ # Computational backends
│ ├── numpy_backend.py # NumPy (CPU)
│ └── cupy_backend.py # CuPy (GPU)
├── io/ # Input/output
│ └── exporters/ # Data export formats
├── utils/ # Utilities
│ └── mode_matching.py # Mode overlap integrals
└── visualization/ # Plotting helpers
Design Patterns
1. Backend Abstraction
All numerical operations go through a backend abstraction:
from prismo.backends import Backend, get_backend
class Solver:
def __init__(self, backend='numpy'):
self.backend = get_backend(backend)
def update_e_fields(self):
# Backend-agnostic code
curl_h = self.backend.curl_h(self.fields.Hx, self.fields.Hy, self.fields.Hz)
self.fields.Ex += self.backend.multiply(self.dt, curl_h.x)
Benefits:
Single codebase for CPU/GPU
Easy to add new backends
Performance portability
2. Component Pattern
Sources, monitors, and boundaries follow a consistent interface:
class Component(ABC):
def initialize(self, grid: YeeGrid) -> None:
"""Initialize component on grid."""
pass
@abstractmethod
def update(self, fields: Fields, time: float, dt: float) -> None:
"""Update component at each time step."""
pass
Example:
class Source(Component):
def update(self, fields, time, dt):
# Add source contribution to fields
pass
class Monitor(Component):
def update(self, fields, time, dt):
# Record field values
pass
3. Yee Grid
The Yee grid staggers E and H fields in space and time:
E fields: (i, j+½, k+½) Ex
(i+½, j, k+½) Ey
(i+½, j+½, k) Ez
H fields: (i+½, j, k) Hx
(i, j+½, k) Hy
(i, j, k+½) Hz
Implementation:
class YeeGrid:
def __init__(self, dimensions, spacing):
self.nx, self.ny, self.nz = dimensions
self.dx, self.dy, self.dz = spacing
def get_component_indices(self, component, region):
"""Get grid indices for a field component."""
# Handle staggering offsets
pass
4. Field Storage
Fields are stored contiguously for cache efficiency:
class ElectromagneticFields:
def __init__(self, grid, backend):
self.Ex = backend.zeros((grid.nx, grid.ny, grid.nz))
self.Ey = backend.zeros((grid.nx, grid.ny, grid.nz))
self.Ez = backend.zeros((grid.nx, grid.ny, grid.nz))
self.Hx = backend.zeros((grid.nx, grid.ny, grid.nz))
self.Hy = backend.zeros((grid.nx, grid.ny, grid.nz))
self.Hz = backend.zeros((grid.nx, grid.ny, grid.nz))
Key Algorithms
FDTD Update Equations
E-field update:
E^{n+1} = E^n + (Δt/ε) ∇ × H^{n+½}
H-field update:
H^{n+½} = H^{n-½} + (Δt/μ) ∇ × E^n
Implementation:
def step(self):
# Update H from E (curl operator)
self.update_h_fields()
# Apply H boundary conditions
self.apply_h_boundaries()
# Update E from H
self.update_e_fields()
# Apply E boundary conditions
self.apply_e_boundaries()
# Add sources
self.apply_sources()
# Record monitors
self.update_monitors()
PML Implementation
Convolutional PML (CPML) for absorbing boundaries:
class CPML:
def __init__(self, thickness, grid):
# Compute PML conductivity profile
self.sigma = self._compute_sigma_profile(thickness)
# Auxiliary fields for convolution
self.Psi_Ex = zeros_like(...)
def update_e_fields(self, Ex, Hy, Hz):
# Update auxiliary fields
self.Psi_Ex = self.b * self.Psi_Ex + self.a * (dHz_dy - dHy_dz)
# Update E with PML correction
Ex += self.dt / self.epsilon * self.Psi_Ex
Mode Solver
Eigenvalue problem for waveguide modes:
def solve_modes(self):
# Build operators
A = self._build_curl_curl_operator() # ∇ × ∇ ×
M = self._build_material_operator() # ε
# Solve: (A + k₀²M) E = β² E
eigenvalues, eigenvectors = sparse_eig(A + k0**2 * M)
# Extract modes
beta = sqrt(eigenvalues)
neff = beta / k0
return modes
Performance Considerations
Memory Layout
Fields are stored in C-contiguous order for cache efficiency:
# Good: Contiguous access
for k in range(nz):
for j in range(ny):
for i in range(nx):
Ex[i, j, k] = ... # Fast
# Bad: Strided access
for i in range(nx):
for j in range(ny):
for k in range(nz):
Ex[i, j, k] = ... # Slower
Vectorization
Use NumPy/CuPy vectorization:
# Good: Vectorized
Ex[1:-1, 1:-1, 1:-1] += dt / eps * (
(Hz[1:-1, 2:, 1:-1] - Hz[1:-1, :-2, 1:-1]) / (2*dy) -
(Hy[1:-1, 1:-1, 2:] - Hy[1:-1, 1:-1, :-2]) / (2*dz)
)
# Bad: Explicit loops
for i in range(1, nx-1):
for j in range(1, ny-1):
for k in range(1, nz-1):
Ex[i,j,k] += ... # Much slower
GPU Considerations
Batch operations: Keep data on GPU, minimize transfers
Kernel fusion: Combine operations to reduce kernel launches
Memory coalescing: Access memory in aligned, contiguous patterns
Extension Points
Adding a New Source
Inherit from
Sourcebase classImplement
update_fields()methodRegister in
sources/__init__.py
from prismo.sources.base import Source
class MySource(Source):
def update_fields(self, fields, time, dt):
# Compute source value
value = self.compute_amplitude(time)
# Add to fields
fields.Ez[self.indices] += value
Adding a New Material Model
Inherit from
MaterialImplement
update_polarization()for dispersive materialsAdd to material library
from prismo.materials.dispersion import DispersiveMaterial
class MyMaterial(DispersiveMaterial):
def update_polarization(self, E, P, time, dt):
# Update polarization based on material model
# E.g., Drude, Lorentz, Debye, etc.
pass
Testing Strategy
# Unit tests: Test individual components
def test_gaussian_pulse():
pulse = GaussianPulse(frequency=1e14, width=1e-15)
assert pulse.value(0.0) == 1.0
assert abs(pulse.value(5e-15)) < 0.01
# Integration tests: Test component interactions
def test_source_monitor_integration():
sim = Simulation(...)
sim.add_source(source)
sim.add_monitor(monitor)
sim.run(...)
assert monitor.has_data()
# Validation tests: Compare to analytical/reference solutions
def test_plane_wave_propagation():
sim = setup_plane_wave_sim()
sim.run(...)
error = compare_to_analytical()
assert error < 0.01
See Also
Contributing to Prismo - How to contribute code
Testing Guide - Testing guidelines
Performance Benchmarks - Performance benchmarks
API Reference - API reference