# Grid1D Complete Tutorial
This tutorial provides comprehensive examples for using the `grid1d` library, from basic usage to advanced techniques.
## Table of Contents
1. [Getting Started](#getting-started)
2. [Understanding Intervals](#understanding-intervals)
3. [Creating Grids](#creating-grids)
4. [Point Location](#point-location)
5. [Working with Different Scalar Types](#working-with-different-scalar-types)
6. [Grid Operations](#grid-operations)
7. [Adaptive Mesh Refinement](#adaptive-mesh-refinement)
8. [Real-World Applications](#real-world-applications)
9. [Performance Optimization](#performance-optimization)
10. [Common Patterns and Best Practices](#common-patterns-and-best-practices)
---
## Getting Started
### Installation
Add to your `Cargo.toml`:
```toml
[dependencies]
grid1d = "0.3.3"
try_create = "0.1.2"
sorted-vec = "0.8"
num-valid = "0.3"
```
### Your First Grid
```rust
use grid1d::{Grid1D, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
fn main() {
// Create a uniform grid from 0 to 10 with 10 equal intervals
let domain = IntervalClosed::new(0.0, 10.0);
let grid = Grid1D::uniform(domain, NumIntervals::try_new(10).unwrap());
println!("Created grid with {} points", grid.coords().len());
// Output: Created grid with 11 points (0, 1, 2, ..., 10)
}
```
**Key takeaways:**
- Domains are defined using interval types (e.g., `IntervalClosed` for `[a, b]`)
- `NumIntervals` represents the number of intervals, not points
- A grid with N intervals has N+1 boundary points
---
## Understanding Intervals
The library provides type-safe interval representations:
```rust
use grid1d::intervals::*;
// Closed interval: [0, 1] - includes both endpoints
let closed = IntervalClosed::new(0.0, 1.0);
assert!(closed.contains_point(&0.0)); // ✓ Included
assert!(closed.contains_point(&1.0)); // ✓ Included
// Open interval: (0, 1) - excludes both endpoints
let open = IntervalOpen::new(0.0, 1.0);
assert!(!open.contains_point(&0.0)); // ✗ Excluded
assert!(!open.contains_point(&1.0)); // ✗ Excluded
// Half-open interval: [0, 1) - includes left, excludes right
let half_open = IntervalLowerClosedUpperOpen::new(0.0, 1.0);
assert!(half_open.contains_point(&0.0)); // ✓ Included
assert!(!half_open.contains_point(&1.0)); // ✗ Excluded
// Unbounded intervals
let unbounded_above = IntervalLowerClosedUpperUnbounded::new(0.0); // [0, +∞)
let unbounded_below = IntervalLowerUnboundedUpperClosed::new(1.0); // (-∞, 1]
```
**Why this matters:**
- Correct mathematical semantics prevent boundary errors
- Type safety ensures you can't confuse open and closed bounds
- Essential for finite element methods and numerical integration
---
## Creating Grids
### Uniform Grids
Perfect for regular discretizations:
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
// Create uniform grid with analytical spacing
let grid = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(100).unwrap()
);
// Properties
assert_eq!(grid.num_intervals().as_ref(), &100);
assert_eq!(grid.coords().len(), 101); // 100 intervals = 101 points
// Spacing is exactly 0.01
let spacing = grid.coords()[1] - grid.coords()[0];
assert!((spacing - 0.01).abs() < 1e-15);
```
### Non-Uniform Grids
For adaptive spacing and complex features:
```rust
use grid1d::{Grid1D, intervals::*};
use sorted_vec::partial::SortedSet;
// Create boundary layer mesh with fine spacing near x=0
let coords = SortedSet::from_unsorted(vec![
0.0, // Wall boundary
0.01, // Fine spacing
0.02,
0.05,
0.1, // Transitional
0.3,
0.6, // Coarse spacing
1.0 // Far boundary
]);
let grid = Grid1D::<IntervalClosed<f64>>::try_from_sorted(
coords
).unwrap();
println!("First interval length: {}",
grid.coords()[1] - grid.coords()[0]); // 0.01
println!("Last interval length: {}",
grid.coords()[7] - grid.coords()[6]); // 0.4
```
### Grid Type Conversions
Convert between specific and generic grid types using `From`/`Into`:
```rust
use grid1d::{Grid1D, Grid1DUniform, Grid1DNonUniform, IntervalPartition, intervals::*, scalars::NumIntervals};
use grid1d::coords::Coords1D;
use sorted_vec::partial::SortedSet;
use try_create::TryNew;
// Create a uniform grid
let uniform_grid = Grid1DUniform::new(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(10).unwrap()
);
// Convert to generic Grid1D using .into()
let generic_grid: Grid1D<IntervalClosed<f64>> = uniform_grid.into();
// Or using From::from()
let uniform_grid2 = Grid1DUniform::new(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(5).unwrap()
);
let generic_grid2 = Grid1D::from(uniform_grid2);
// Same works for non-uniform grids
let coords = Coords1D::try_from(
SortedSet::from_unsorted(vec![0.0, 0.3, 0.7, 1.0])
).unwrap();
let non_uniform = Grid1DNonUniform::try_new(coords).unwrap();
let generic_from_non_uniform: Grid1D<_> = non_uniform.into();
```
**Use cases for conversions:**
- **Generic functions**: Write functions that accept `Grid1D<D>` to work with any grid type
- **Collections**: Store different grid types in the same `Vec<Grid1D<D>>`
- **API boundaries**: Return generic `Grid1D` from functions that may create either type internally
### Using Utility Functions
```rust
use grid1d::{intervals::*, scalars::NumIntervals, linspace, logspace};
use sorted_vec::partial::SortedSet;
use try_create::TryNew;
// Linear spacing
let domain = IntervalClosed::new(0.0, 1.0);
let n = NumIntervals::try_new(10).unwrap();
let linear_coords = linspace(&domain, &n);
// Logarithmic spacing (useful for wide-range problems)
// logspace takes: base, exponent_domain, num_intervals
// For points from 10^0=1 to 10^3=1000:
let exponent_domain = IntervalClosed::new(0.0, 3.0); // exponents: 0, 1, 2, 3
let log_coords = logspace(&10.0_f64, &exponent_domain, &NumIntervals::try_new(3).unwrap());
// Result: [1.0, 10.0, 100.0, 1000.0]
// Convert to grid (domain inferred from coordinates)
let log_grid = grid1d::Grid1D::<IntervalClosed<f64>>::try_from_sorted(
SortedSet::from_unsorted(log_coords)
).unwrap();
```
---
## Point Location
One of the core operations is finding which interval contains a given point.
### Basic Point Location
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
let grid = Grid1D::uniform(
IntervalClosed::new(0.0, 10.0),
NumIntervals::try_new(10).unwrap()
);
// Find interval containing a point
let point = 3.7;
let interval_id = grid.find_interval_id_of_point(&point);
println!("Point {} is in interval {}", point, interval_id.as_ref());
// Access the interval
let interval = grid.interval(&interval_id);
println!("Interval: {:?}", interval);
// Get interval properties
let length = grid.interval_length(&interval_id);
println!("Interval length: {}", length.as_ref());
```
### Batch Point Location
For efficiency when locating many points:
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
let grid = Grid1D::uniform(
IntervalClosed::new(0.0, 10.0),
NumIntervals::try_new(100).unwrap()
);
// Locate multiple points at once
let points = vec![1.5, 3.7, 5.2, 7.8, 9.1];
let intervals = grid.find_intervals_for_points(&points);
for (point, interval_id) in points.iter().zip(intervals.iter()) {
println!("Point {} → Interval {}", point, interval_id.as_ref());
}
```
### Parallel Batch Point Location
For large-scale point location, use the parallel version with Rayon:
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
let grid = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(10_000).unwrap()
);
// Generate many points to locate
let points: Vec<f64> = (0..100_000)
.map(|i| (i as f64) / 100_000.0)
.collect();
// Parallel point location using Rayon
let intervals = grid.find_intervals_for_points_parallel(&points);
assert_eq!(intervals.len(), points.len());
```
**When to use parallel vs sequential:**
| < 1,000 | Sequential (`find_intervals_for_points`) |
| ≥ 1,000 | Parallel (`find_intervals_for_points_parallel`) |
> **Note**: The parallel version requires the grid type to be `Sync` and the point type to be `Send + Sync`.
> All built-in grid types satisfy these requirements.
### Performance Comparison
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use sorted_vec::partial::SortedSet;
use try_create::TryNew;
use std::time::Instant;
// Uniform grid: O(1) point location
let uniform_grid = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(100_000).unwrap()
);
let start = Instant::now();
for i in 0..10_000 {
let point = i as f64 / 10_000.0;
let _ = uniform_grid.find_interval_id_of_point(&point);
}
println!("Uniform grid: {:?}", start.elapsed());
// Non-uniform grid: O(log n) point location
let coords: Vec<f64> = (0..=100_000)
.map(|i| (i as f64 / 100_000.0).powi(2)) // Non-uniform spacing
.collect();
let non_uniform_grid = Grid1D::<IntervalClosed<f64>>::try_from_sorted(
SortedSet::from_unsorted(coords)
).unwrap();
let start = Instant::now();
for i in 0..10_000 {
let point = (i as f64 / 10_000.0).powi(2);
let _ = non_uniform_grid.find_interval_id_of_point(&point);
}
println!("Non-uniform grid: {:?}", start.elapsed());
```
---
## Working with Different Scalar Types
### Choosing the Right Scalar Type
```rust
use grid1d::{Grid1D, intervals::*, scalars::NumIntervals};
use num_valid::{RealNative64StrictFiniteInDebug, RealNative64StrictFinite};
use try_create::TryNew;
// Option 1: Raw f64 (maximum performance, no safety checks)
let fast_grid = Grid1D::uniform(
IntervalClosed::new(0.0_f64, 1.0_f64),
NumIntervals::try_new(100).unwrap()
);
// Option 2: Debug-validated (recommended - f64 performance in release)
type Real = RealNative64StrictFiniteInDebug;
let safe_grid = Grid1D::uniform(
IntervalClosed::new(
Real::try_new(0.0).unwrap(),
Real::try_new(1.0).unwrap()
),
NumIntervals::try_new(100).unwrap()
);
// Option 3: Always validated (safety-critical applications)
type SafeReal = RealNative64StrictFinite;
let validated_grid = Grid1D::uniform(
IntervalClosed::new(
SafeReal::try_new(0.0).unwrap(),
SafeReal::try_new(1.0).unwrap()
),
NumIntervals::try_new(100).unwrap()
);
```
### Arbitrary Precision Arithmetic
For high-precision scientific computing:
```rust
#[cfg(feature = "rug")]
use grid1d::{Grid1D, intervals::*, scalars::NumIntervals};
#[cfg(feature = "rug")]
use num_valid::RealRugStrictFinite;
#[cfg(feature = "rug")]
use try_create::TryNew;
#[cfg(feature = "rug")]
fn high_precision_example() {
// 256-bit precision
type HighPrecReal = RealRugStrictFinite<256>;
let hp_grid = Grid1D::uniform(
IntervalClosed::new(
HighPrecReal::try_new(0.0).unwrap(),
HighPrecReal::try_new(1.0).unwrap()
),
NumIntervals::try_new(100).unwrap()
);
// All operations work identically to f64 version
let hp_point = HighPrecReal::try_new(0.333333333333333333333).unwrap();
let interval_id = hp_grid.find_interval_id_of_point(&hp_point);
}
```
---
## Grid Operations
### Grid Intersection
Find intervals that overlap with a subdomain:
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
let grid = Grid1D::uniform(
IntervalClosed::new(0.0, 10.0),
NumIntervals::try_new(10).unwrap()
);
// Find intervals intersecting subdomain [3, 7]
let subdomain = IntervalClosed::new(3.0, 7.0);
let intersecting = grid.intervals_in_intersection(&subdomain);
println!("Intervals intersecting [{}, {}]:",
subdomain.lower_bound_value(),
subdomain.upper_bound_value());
for interval_id in intersecting {
println!(" Interval {}", interval_id.as_ref());
}
```
### Grid Union for Multi-Physics
Combine grids from different physics or solvers:
```rust
use grid1d::{Grid1D, Grid1DUnion, intervals::*, scalars::NumIntervals};
use sorted_vec::partial::SortedSet;
use try_create::TryNew;
let domain = IntervalClosed::new(0.0, 10.0);
// Coarse global grid for fluid flow
let fluid_grid = Grid1D::uniform(domain.clone(), NumIntervals::try_new(10).unwrap());
// Fine local grid for chemistry
let chemistry_coords = SortedSet::from_unsorted(vec![
0.0, 2.0, 2.1, 2.2, 2.5, 3.0, 3.1, 3.2, 5.0, 10.0
]);
let chemistry_grid = Grid1D::<IntervalClosed<f64>>::try_from_sorted(chemistry_coords).unwrap();
// Create unified grid
let unified = Grid1DUnion::try_new(&fluid_grid, &chemistry_grid).unwrap();
println!("Fluid grid: {} intervals", fluid_grid.num_intervals().as_ref());
println!("Chemistry grid: {} intervals", chemistry_grid.num_intervals().as_ref());
println!("Unified grid: {} intervals", unified.num_refined_intervals().as_ref());
// Map data between grids
for (unified_id, fluid_id, chem_id) in unified.iter_interval_mappings() {
println!("Unified[{}] ← Fluid[{}], Chemistry[{}]",
unified_id.as_ref(), fluid_id.as_ref(), chem_id.as_ref());
}
```
---
## Adaptive Mesh Refinement
### Uniform Refinement
Double the resolution everywhere:
```rust
use grid1d::{Grid1D, intervals::*, scalars::{NumIntervals, PositiveNumPoints1D}};
use try_create::TryNew;
let base_grid = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(4).unwrap()
);
// Add 1 point per interval (doubles resolution)
let refinement = base_grid.refine_uniform(&PositiveNumPoints1D::try_new(1).unwrap());
println!("Base grid: {} intervals", base_grid.num_intervals().as_ref());
println!("Refined grid: {} intervals",
refinement.refined_grid().num_intervals().as_ref());
// Access refined grid
let refined = refinement.refined_grid();
```
### Selective Refinement
Refine only where needed (adaptive mesh refinement):
```rust
use grid1d::{Grid1D, intervals::*, scalars::*};
use try_create::TryNew;
let base_grid = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(10).unwrap()
);
// Define refinement plan: (interval_id, num_additional_points)
use std::collections::BTreeMap;
let refinement_plan = BTreeMap::from([
(IntervalId::new(2), PositiveNumPoints1D::try_new(3).unwrap()), // 4 sub-intervals
(IntervalId::new(3), PositiveNumPoints1D::try_new(3).unwrap()), // 4 sub-intervals
(IntervalId::new(7), PositiveNumPoints1D::try_new(1).unwrap()), // 2 sub-intervals
]);
let refined = base_grid.refine(&refinement_plan);
println!("Refined grid has {} intervals",
refined.refined_grid().num_intervals().as_ref());
// Trace refinement lineage
for (refined_id, base_id) in refined.iter_refined_with_mapping() {
println!("Refined interval {} came from base interval {}",
refined_id.as_ref(), base_id.as_ref());
}
```
---
## Real-World Applications
### Finite Difference Methods
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
use std::ops::Deref;
fn solve_heat_equation() {
// Setup spatial grid
let grid = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(100).unwrap()
);
let coords = grid.coords().deref();
let n = coords.len();
let dx = coords[1] - coords[0];
let dt = 0.0001; // Time step
let alpha = 1.0; // Thermal diffusivity
// Stability criterion: dt ≤ dx²/(2α)
assert!(dt <= dx * dx / (2.0 * alpha));
// Initialize solution
let mut u = vec![0.0; n];
// Initial condition: u(x, 0) = sin(πx)
for (i, &x) in coords.iter().enumerate() {
u[i] = (std::f64::consts::PI * x).sin();
}
// Time stepping
let num_steps = 1000;
let mut u_new = vec![0.0; n];
for _step in 0..num_steps {
// Boundary conditions: u(0,t) = u(1,t) = 0
u_new[0] = 0.0;
u_new[n-1] = 0.0;
// Interior points: forward Euler
for i in 1..n-1 {
u_new[i] = u[i] + alpha * dt / (dx * dx) *
(u[i-1] - 2.0*u[i] + u[i+1]);
}
// Update
std::mem::swap(&mut u, &mut u_new);
}
println!("Final temperature at x=0.5: {}", u[n/2]);
}
```
### Boundary Layer Mesh Generation
```rust
use grid1d::{Grid1D, intervals::*};
use sorted_vec::partial::SortedSet;
fn create_boundary_layer_mesh(
length: f64,
boundary_thickness: f64,
n_boundary: usize,
n_outer: usize
) -> Grid1D<IntervalClosed<f64>> {
let mut points = vec![0.0];
// Boundary layer: quadratic clustering
for i in 1..=n_boundary {
let eta = (i as f64) / (n_boundary as f64);
let y = boundary_thickness * eta * eta;
points.push(y);
}
// Outer region: linear spacing
for i in 1..=n_outer {
let y = boundary_thickness +
(length - boundary_thickness) * (i as f64) / (n_outer as f64);
points.push(y);
}
Grid1D::<IntervalClosed<f64>>::try_from_sorted(
SortedSet::from_unsorted(points)
).unwrap()
}
fn main() {
let cfd_grid = create_boundary_layer_mesh(1.0, 0.05, 30, 20);
println!("CFD grid: {} intervals", cfd_grid.num_intervals().as_ref());
// First interval is very fine
let coords = cfd_grid.coords();
println!("First interval: Δy = {:.6}", coords[1] - coords[0]);
// Last interval is coarse
let n = coords.len();
println!("Last interval: Δy = {:.6}", coords[n-1] - coords[n-2]);
}
```
---
## Performance Optimization
### Choosing Grid Type for Performance
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
use std::time::Instant;
fn benchmark_grid_types() {
let n_intervals = 100_000;
let n_queries = 10_000;
// Uniform grid: O(1) lookup
let uniform = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(n_intervals).unwrap()
);
let start = Instant::now();
for i in 0..n_queries {
let point = (i as f64) / (n_queries as f64);
let _ = uniform.find_interval_id_of_point(&point);
}
let uniform_time = start.elapsed();
println!("Uniform grid (O(1)): {:?}", uniform_time);
println!("Average per query: {:?}", uniform_time / n_queries);
}
```
### Batch Operations
```rust
use grid1d::{Grid1D, IntervalPartition, intervals::*, scalars::NumIntervals};
use try_create::TryNew;
fn batch_point_location() {
let grid = Grid1D::uniform(
IntervalClosed::new(0.0, 1.0),
NumIntervals::try_new(1000).unwrap()
);
// Generate many query points
let points: Vec<f64> = (0..10000)
.map(|i| (i as f64) / 10000.0)
.collect();
// Batch locate (more efficient than loop)
let intervals = grid.find_intervals_for_points(&points);
assert_eq!(points.len(), intervals.len());
}
```
---
## Common Patterns and Best Practices
### Pattern 1: Type Alias for Scalar Types
```rust
use grid1d::{Grid1D, intervals::*, scalars::NumIntervals};
use num_valid::RealNative64StrictFiniteInDebug;
use try_create::TryNew;
// Define once at module level
type Real = RealNative64StrictFiniteInDebug;
type Domain = IntervalClosed<Real>;
fn create_grid(a: f64, b: f64, n: usize) -> Grid1D<Domain> {
Grid1D::uniform(
IntervalClosed::new(
Real::try_new(a).unwrap(),
Real::try_new(b).unwrap()
),
NumIntervals::try_new(n).unwrap()
)
}
```
### Pattern 2: Error Handling
```rust
use grid1d::{Grid1D, intervals::*, ErrorsGrid1D};
use sorted_vec::partial::SortedSet;
fn create_grid_from_user_input(coords: Vec<f64>) -> Result<Grid1D<IntervalClosed<f64>>, ErrorsGrid1D<IntervalClosed<f64>>> {
let sorted_coords = SortedSet::from_unsorted(coords);
// Domain is inferred from the first and last coordinates
Grid1D::try_from_sorted(sorted_coords)
}
```
### Pattern 3: Generic Functions Over Grid Types
```rust
use grid1d::{IntervalPartition, intervals::IntervalTrait};
use num_valid::RealScalar;
fn analyze_grid<G>(grid: &G) -> String
where
G: IntervalPartition,
G::Domain1D: IntervalTrait,
<G::Domain1D as IntervalTrait>::RealType: RealScalar + std::fmt::Display,
{
format!(
"Grid with {} intervals spanning domain {:?}",
grid.num_intervals().as_ref(),
grid.domain()
)
}
```
### Pattern 4: Reusable Grid Generators
```rust
use grid1d::{Grid1D, intervals::*, scalars::NumIntervals, linspace};
use sorted_vec::partial::SortedSet;
use try_create::TryNew;
/// Generates grids with various standard distributions
enum GridDistribution {
Uniform,
Logarithmic,
Geometric { ratio: f64 },
Custom { coords: Vec<f64> },
}
fn generate_grid(
domain: IntervalClosed<f64>,
n: usize,
dist: GridDistribution
) -> Grid1D<IntervalClosed<f64>> {
let num_intervals = NumIntervals::try_new(n).unwrap();
match dist {
GridDistribution::Uniform => {
Grid1D::uniform(domain, num_intervals)
}
GridDistribution::Custom { coords } => {
Grid1D::<IntervalClosed<f64>>::try_from_sorted(SortedSet::from_unsorted(coords)).unwrap()
}
// Implement other distributions...
_ => Grid1D::uniform(domain, num_intervals)
}
}
```
---
## Conclusion
This tutorial covers the essential usage patterns of the `grid1d` library. For more advanced topics:
- **API Documentation**: <https://docs.rs/grid1d>
- **Source Code**: <https://gitlab.com/max.martinelli/grid1d>
- **Examples**: See the `examples/` directory
### Next Steps
1. **Start simple**: Use uniform grids with `f64` for learning
2. **Add safety**: Switch to `RealNative64StrictFiniteInDebug` for production
3. **Optimize**: Use specialized grid types based on profiling
4. **Scale up**: Try arbitrary precision for high-accuracy requirements
Happy gridding! 🎯