diff --git a/crates/clawhdf5-format/src/chunked_read.rs b/crates/clawhdf5-format/src/chunked_read.rs index 70176e2..7d3efd0 100644 --- a/crates/clawhdf5-format/src/chunked_read.rs +++ b/crates/clawhdf5-format/src/chunked_read.rs @@ -6,15 +6,16 @@ extern crate alloc; #[cfg(not(feature = "std"))] use alloc::{format, vec, vec::Vec}; -use crate::chunk_cache::CacheAlignedBuffer; #[cfg(feature = "std")] -use crate::chunk_cache::ChunkCache; +use crate::chunk_cache::{CacheAlignedBuffer, ChunkCache}; use crate::data_layout::DataLayout; use crate::dataspace::Dataspace; use crate::datatype::Datatype; use crate::error::FormatError; use crate::extensible_array::{ExtensibleArrayHeader, read_extensible_array_chunks}; use crate::filter_pipeline::FilterPipeline; +use crate::filters::{DecodeScratch, decompress_chunk_exact_with}; +#[cfg(feature = "std")] use crate::filters::{all_filters_skipped, decompress_chunk_exact}; use crate::fixed_array::{FixedArrayHeader, read_fixed_array_chunks}; #[cfg(feature = "std")] @@ -26,60 +27,329 @@ use crate::parallel_read; #[cfg(feature = "parallel")] use crate::lane_partition::PartitionStats; -/// Decompress all chunks into cache-line-aligned buffers, using lane-partitioned -/// parallel decompression when the `parallel` feature is enabled and the chunk -/// count exceeds the threshold. -fn decompress_all_chunks( - file_data: &[u8], - chunks: &[ChunkInfo], - pipeline: Option<&FilterPipeline>, - chunk_total_bytes: usize, - element_size: u32, -) -> Result, FormatError> { - #[cfg(feature = "parallel")] +/// Run `f` with this thread's chunk-decoding scratch buffers (see +/// [`DecodeScratch`]), kept between reads so decoding reuses memory instead +/// of faulting in fresh pages for every chunk. A re-entrant call (a +/// registered filter codec that itself reads a file) gets a fresh scratch. +fn with_scratch(f: impl FnOnce(&mut DecodeScratch) -> R) -> R { + #[cfg(feature = "std")] { - if let Some(pl) = pipeline - && parallel_read::should_use_parallel(chunks.len()) - && parallel_read::pool_can_parallelise() - { - // Seed from the first chunk's address and count for determinism. - let seed = chunks.first().map(|c| c.address).unwrap_or(0) ^ (chunks.len() as u64); - let (data, _stats) = parallel_read::decompress_chunks_lane_partitioned( - file_data, - chunks, - pl, - chunk_total_bytes, - element_size, - seed, - None, // auto-detect lane count - )?; - return Ok(data.into_iter().map(CacheAlignedBuffer::from_vec).collect()); + use std::cell::RefCell; + std::thread_local! { + static SCRATCH: RefCell = RefCell::new(DecodeScratch::new()); + } + let mut f = Some(f); + let kept = SCRATCH.try_with(|cell| { + let mut scratch = cell.try_borrow_mut().ok()?; + let f = f.take()?; + let r = f(&mut scratch); + scratch.trim(); + Some(r) + }); + if let Ok(Some(r)) = kept { + return r; + } + let f = f.expect("with_scratch: closure already run"); + f(&mut DecodeScratch::new()) + } + #[cfg(not(feature = "std"))] + f(&mut DecodeScratch::new()) +} + +/// A full read's output buffer, written through a raw pointer so that +/// several threads can place chunks into it at once. +struct OutBuf<'a> { + ptr: *mut u8, + len: usize, + _borrow: core::marker::PhantomData<&'a mut [u8]>, +} + +// SAFETY: `OutBuf` is a `&mut [u8]` that hands out writes; sharing it across +// threads is sound as long as concurrent writes do not overlap, which +// `write`'s contract requires. +unsafe impl Send for OutBuf<'_> {} +unsafe impl Sync for OutBuf<'_> {} + +impl<'a> OutBuf<'a> { + fn new(out: &'a mut [u8]) -> Self { + Self { + ptr: out.as_mut_ptr(), + len: out.len(), + _borrow: core::marker::PhantomData, } } - // Sequential fallback — allocate into aligned buffers - let mut result = Vec::with_capacity(chunks.len()); - for chunk_info in chunks { - let c_addr = chunk_info.address as usize; - let size = chunk_info.chunk_size as usize; - ensure_len(file_data, c_addr, size)?; - let raw_chunk = &file_data[c_addr..c_addr + size]; - - let decompressed = if let Some(pl) = pipeline { - decompress_chunk_exact( - raw_chunk, - pl, - chunk_total_bytes, - element_size, - chunk_info.filter_mask, - &chunk_info.offsets, - )? - } else { - raw_chunk.to_vec() - }; - result.push(CacheAlignedBuffer::from_vec(decompressed)); + fn len(&self) -> usize { + self.len } - Ok(result) + + /// Copy `src` to `[at, at + src.len())`. Out of range is a no-op (the + /// callers check first). + /// + /// # Safety + /// + /// No other thread may be writing an overlapping range at the same time. + unsafe fn write(&self, at: usize, src: &[u8]) { + if at.checked_add(src.len()).is_none_or(|end| end > self.len) { + debug_assert!(false, "OutBuf::write out of range"); + return; + } + // SAFETY: in range (checked above) of a live `&'a mut [u8]`; `src` + // cannot overlap it (the output is exclusively borrowed); no + // concurrent overlapping write (the caller's contract). + unsafe { core::ptr::copy_nonoverlapping(src.as_ptr(), self.ptr.add(at), src.len()) } + } +} + +/// How a chunked dataset's chunks map into its row-major output. +struct ChunkPlacer { + rank: usize, + chunk_dims: Vec, + ds_dims: Vec, + ds_strides: Vec, + chunk_strides: Vec, + elem_size: usize, +} + +impl ChunkPlacer { + /// `chunk_dims` and `ds_dims` have one entry per dimension; the dataset + /// must not be empty (so the stride products stay in range). + fn new(chunk_dims: &[usize], ds_dims: &[usize], elem_size: usize) -> Self { + let rank = chunk_dims.len(); + let mut ds_strides = vec![1usize; rank]; + let mut chunk_strides = vec![1usize; rank]; + for i in (0..rank.saturating_sub(1)).rev() { + ds_strides[i] = ds_strides[i + 1] * ds_dims[i + 1]; + chunk_strides[i] = chunk_strides[i + 1] * chunk_dims[i + 1]; + } + Self { + rank, + chunk_dims: chunk_dims.to_vec(), + ds_dims: ds_dims.to_vec(), + ds_strides, + chunk_strides, + elem_size, + } + } + + /// Copy the decoded chunk `data`, whose first element is at `offsets` + /// (one per dimension), to its place in `out`. + /// + /// # Safety + /// + /// No other thread may be placing a chunk whose region overlaps this + /// one's (see [`Self::regions_disjoint`]). + unsafe fn place(&self, data: &[u8], offsets: &[u64], out: &OutBuf<'_>) { + if self.rank == 0 { + let copy_len = data.len().min(out.len()); + // SAFETY: in range; exclusivity is the caller's contract. + unsafe { out.write(0, &data[..copy_len]) }; + return; + } + let mut short = [0usize; 8]; + let mut long = Vec::new(); + let chunk_offsets: &mut [usize] = if self.rank <= short.len() { + &mut short[..self.rank] + } else { + long.resize(self.rank, 0); + &mut long + }; + for (o, &off) in chunk_offsets.iter_mut().zip(offsets) { + *o = off as usize; + } + // SAFETY: the caller's contract. + unsafe { + copy_chunk_into( + data, + out, + chunk_offsets, + &self.chunk_dims, + &self.ds_dims, + &self.ds_strides, + &self.chunk_strides, + self.elem_size, + self.rank, + ) + }; + } + + /// Whether placing `chunks` writes pairwise disjoint regions of the + /// output, so they can be placed concurrently: every chunk starts on the + /// chunk grid and no two start at the same place (a chunk wholly outside + /// the dataset writes nothing and is ignored). A corrupt index that + /// breaks this is read one chunk at a time instead. + #[cfg(feature = "parallel")] + fn regions_disjoint(&self, chunks: &[ChunkInfo]) -> bool { + if self.rank == 0 { + return chunks.len() <= 1; + } + let mut keys = Vec::with_capacity(chunks.len()); + 'chunks: for c in chunks { + if c.offsets.len() < self.rank { + return false; + } + let mut key = 0u64; + for d in 0..self.rank { + let (off, cd, dd) = ( + c.offsets[d], + self.chunk_dims[d] as u64, + self.ds_dims[d] as u64, + ); + if off >= dd { + continue 'chunks; + } + if off % cd != 0 { + return false; + } + let Some(k) = key + .checked_mul(dd.div_ceil(cd)) + .and_then(|k| k.checked_add(off / cd)) + else { + return false; + }; + key = k; + } + keys.push(key); + } + keys.sort_unstable(); + keys.windows(2).all(|w| w[0] != w[1]) + } +} + +/// The per-file chunk cache a full read uses, if any: the cache, this +/// dataset's key in it (its chunk-index address), and whether the dataset +/// fits in it — only then are its chunks looked up and inserted. +#[cfg(feature = "std")] +type CacheUse<'a> = Option<(&'a ChunkCache, u64, bool)>; +#[cfg(not(feature = "std"))] +type CacheUse<'a> = Option<&'a core::convert::Infallible>; + +/// Decode every chunk in `chunks` and place it in `output` (the dataset's +/// whole row-major extent, zeroed). +/// +/// Each chunk goes straight from its decoder to its place in the output: +/// decoded into this thread's reusable scratch (or, when the cache is to +/// keep it, into a buffer the cache takes), then copied. With the +/// `parallel` feature and a filter pipeline, the chunks are shared out +/// between the calling thread and idle rayon workers +/// ([`parallel_read::run_with_helpers`]), so the caller never waits on a +/// busy pool. The error returned is the first failing chunk's, in `chunks` +/// order. +#[allow(clippy::too_many_arguments)] +fn fill_from_chunks( + file_data: &[u8], + chunks: &[ChunkInfo], + pipeline: Option<&FilterPipeline>, + placer: &ChunkPlacer, + chunk_total_bytes: usize, + cache: CacheUse<'_>, + output: &mut [u8], +) -> Result<(), FormatError> { + let rank = placer.rank; + let elem_size = placer.elem_size as u32; + let out = OutBuf::new(output); + #[cfg(not(feature = "std"))] + let _ = cache; + + // Decode chunk `i` and place it. Its callers below run it either on one + // thread, or on several for chunks whose regions are pairwise disjoint, + // each chunk once: no two threads ever write the same bytes. + let work = |i: usize, scratch: &mut DecodeScratch| -> Result<(), FormatError> { + let c = &chunks[i]; + if c.offsets.len() < rank { + return Err(FormatError::ChunkedReadError(format!( + "chunk index entry has {} offsets for a rank-{rank} dataset", + c.offsets.len() + ))); + } + let offsets = &c.offsets[..rank]; + let c_addr = c.address as usize; + let size = c.chunk_size as usize; + ensure_len(file_data, c_addr, size)?; + let raw = &file_data[c_addr..c_addr + size]; + let Some(pl) = pipeline else { + // SAFETY: see above. + unsafe { placer.place(raw, offsets, &out) }; + return Ok(()); + }; + // A chunk stored as-is (every filter skipped) is checked and placed + // straight from the file bytes, never cached. + #[cfg(feature = "std")] + if let Some((cache, key, true)) = cache + && !all_filters_skipped(pl, c.filter_mask) + { + let cached = match cache.get_decompressed_in(key, offsets) { + Some(hit) => hit, + None => { + let data = decompress_chunk_exact( + raw, + pl, + chunk_total_bytes, + elem_size, + c.filter_mask, + &c.offsets, + )?; + cache.put_decompressed_in(key, offsets.to_vec(), data) + } + }; + // SAFETY: see above. + unsafe { placer.place(&cached, offsets, &out) }; + return Ok(()); + } + let data = decompress_chunk_exact_with( + raw, + pl, + chunk_total_bytes, + elem_size, + c.filter_mask, + &c.offsets, + scratch, + )?; + // SAFETY: see above. + unsafe { placer.place(data, offsets, &out) }; + Ok(()) + }; + + #[cfg(feature = "parallel")] + if pipeline.is_some() + && parallel_read::should_use_parallel(chunks.len()) + && placer.regions_disjoint(chunks) + { + use core::sync::atomic::{AtomicUsize, Ordering}; + use std::sync::{Mutex, PoisonError}; + + let n = chunks.len(); + let next = AtomicUsize::new(0); + let failed: Mutex> = Mutex::new(None); + let body = || { + with_scratch(|scratch| { + loop { + let i = next.fetch_add(1, Ordering::Relaxed); + if i >= n { + break; + } + if let Err(e) = work(i, scratch) { + let mut failed = failed.lock().unwrap_or_else(PoisonError::into_inner); + if failed.as_ref().is_none_or(|(at, _)| i < *at) { + *failed = Some((i, e)); + } + // Stop handing out chunks. Every chunk before `i` was + // claimed already and finishes, so the error kept is + // the first in order. + next.store(n, Ordering::Relaxed); + break; + } + } + }) + }; + parallel_read::run_with_helpers(parallel_read::helper_count(n), &body); + return match failed.into_inner().unwrap_or_else(PoisonError::into_inner) { + Some((_, e)) => Err(e), + None => Ok(()), + }; + } + + with_scratch(|scratch| (0..chunks.len()).try_for_each(|i| work(i, scratch))) } /// Decompress all chunks with lane-partitioned parallelism and return @@ -557,11 +827,6 @@ pub fn generate_implicit_chunks( chunks } -/// Read a chunked dataset, decompressing chunks as needed. -/// Chunks decompressed together before being copied out, bounding the extra -/// memory a parallel full read holds at once. -const DECODE_BATCH: usize = 128; - /// B-tree v2 record types used for chunk indexing. const BT2_CHUNK_UNFILTERED: u8 = 10; const BT2_CHUNK_FILTERED: u8 = 11; @@ -820,6 +1085,111 @@ pub fn list_chunks( Ok((chunks, chunk_dims)) } +/// The chunk cache a full read may use (`None` without `std`). +#[cfg(feature = "std")] +pub(crate) type CacheRef<'a> = Option<&'a ChunkCache>; +/// The chunk cache a full read may use (`None` without `std`). +#[cfg(not(feature = "std"))] +pub(crate) type CacheRef<'a> = Option<&'a core::convert::Infallible>; + +/// The body of every full chunked read: list the chunks (through the +/// cache's index for this dataset when there is a cache), allocate the +/// output with `alloc` (zeroed, `total_bytes` long, as bytes through +/// `bytes`), and decode every chunk straight into it. +#[allow(clippy::too_many_arguments)] +pub(crate) fn read_chunked_full( + file_data: &[u8], + layout: &DataLayout, + dataspace: &Dataspace, + datatype: &Datatype, + pipeline: Option<&FilterPipeline>, + offset_size: u8, + length_size: u8, + cache: CacheRef<'_>, + alloc: impl FnOnce(usize) -> Result, + bytes: impl FnOnce(&mut O) -> &mut [u8], +) -> Result { + check_chunk_element_size(layout, datatype, offset_size)?; + let elem_size = datatype.type_size() as usize; + let list = || { + list_chunks( + file_data, + layout, + dataspace, + elem_size, + offset_size, + length_size, + ) + }; + + #[cfg(feature = "std")] + let (chunks, chunk_dims, cache_key) = match cache { + Some(cache) => { + let (chunk_dimensions, version, addr_opt) = match layout { + DataLayout::Chunked { + chunk_dimensions, + version, + btree_address, + .. + } => (chunk_dimensions, *version, *btree_address), + _ => { + return Err(FormatError::ChunkedReadError( + "expected chunked layout".into(), + )); + } + }; + let addr = addr_opt.ok_or_else(|| { + FormatError::ChunkedReadError("no address for chunked layout".into()) + })?; + let (rank, chunk_dims) = + chunk_geometry(chunk_dimensions, version, dataspace, elem_size)?; + // The per-file cache is shared across datasets (and threads); + // every lookup is keyed by this dataset's chunk-index address, so + // another dataset's index or chunks are never used for this read. + let chunks = cache.chunks_for(addr, rank, || list().map(|(chunks, _)| chunks))?; + (chunks, chunk_dims, Some((cache, addr))) + } + None => { + let (chunks, chunk_dims) = list()?; + (chunks, chunk_dims, None) + } + }; + #[cfg(not(feature = "std"))] + let (chunks, chunk_dims) = { + let _ = cache; + list()? + }; + + let total_bytes = checked_byte_len(dataspace.checked_num_elements()?, elem_size)?; + let mut output = alloc(total_bytes)?; + if total_bytes == 0 { + // Also keeps the stride products in range: with a zero-sized + // dimension the total is 0 even if other dimensions are huge. + return Ok(output); + } + let ds_dims: Vec = dataspace.dimensions.iter().map(|&d| d as usize).collect(); + let placer = ChunkPlacer::new(&chunk_dims, &ds_dims, elem_size); + let chunk_total_bytes = checked_chunk_byte_len(&chunk_dims, elem_size)?; + // Chunks are cached only when the whole dataset fits: pushing a larger + // dataset through the cache just evicts each chunk moments after + // inserting it. + #[cfg(feature = "std")] + let cache_use = cache_key.map(|(cache, key)| (cache, key, total_bytes <= cache.max_bytes())); + #[cfg(not(feature = "std"))] + let cache_use = None; + fill_from_chunks( + file_data, + &chunks, + pipeline, + &placer, + chunk_total_bytes, + cache_use, + bytes(&mut output), + )?; + Ok(output) +} + +/// Read a chunked dataset, decompressing chunks as needed. pub fn read_chunked_data( file_data: &[u8], layout: &DataLayout, @@ -829,118 +1199,26 @@ pub fn read_chunked_data( offset_size: u8, length_size: u8, ) -> Result, FormatError> { - check_chunk_element_size(layout, datatype, offset_size)?; - let elem_size = datatype.type_size() as usize; - let (chunks, chunk_dims) = list_chunks( + read_chunked_full( file_data, layout, dataspace, - elem_size, + datatype, + pipeline, offset_size, length_size, - )?; - let rank = chunk_dims.len(); - let ds_dims: Vec = dataspace.dimensions.iter().map(|&d| d as usize).collect(); - - // Assemble output - let total_bytes = checked_byte_len(dataspace.checked_num_elements()?, elem_size)?; - if total_bytes == 0 { - // Also keeps the stride products below in range: with a zero-sized - // dimension the total is 0 even if other dimensions are huge. - return Ok(Vec::new()); - } - let mut output = alloc_output(total_bytes)?; - - let mut ds_strides = vec![1usize; rank]; - for i in (0..rank.saturating_sub(1)).rev() { - ds_strides[i] = ds_strides[i + 1] * ds_dims[i + 1]; - } - - let mut chunk_strides = vec![1usize; rank]; - for i in (0..rank.saturating_sub(1)).rev() { - chunk_strides[i] = chunk_strides[i + 1] * chunk_dims[i + 1]; - } - - let chunk_total_bytes = checked_chunk_byte_len(&chunk_dims, elem_size)?; - - // Fast path: no filters — copy directly from file_data without intermediate alloc - if pipeline.is_none() { - for chunk_info in &chunks { - let chunk_offsets: Vec = chunk_info - .offsets - .iter() - .take(rank) - .map(|&o| o as usize) - .collect(); - - let c_addr = chunk_info.address as usize; - let size = chunk_info.chunk_size as usize; - ensure_len(file_data, c_addr, size)?; - let chunk_data = &file_data[c_addr..c_addr + size]; - - if rank == 0 { - let copy_len = chunk_data.len().min(output.len()); - output[..copy_len].copy_from_slice(&chunk_data[..copy_len]); - } else { - copy_chunk_to_output( - chunk_data, - &mut output, - &chunk_offsets, - &chunk_dims, - &ds_dims, - &ds_strides, - &chunk_strides, - elem_size, - rank, - ); - } - } - return Ok(output); - } - - // Filtered path: decompress all chunks then assemble - let decompressed_chunks = decompress_all_chunks( - file_data, - &chunks, - pipeline, - chunk_total_bytes, - elem_size as u32, - )?; - - for (chunk_info, decompressed) in chunks.iter().zip(decompressed_chunks.iter()) { - let chunk_offsets: Vec = chunk_info - .offsets - .iter() - .take(rank) - .map(|&o| o as usize) - .collect(); - - if rank == 0 { - let copy_len = decompressed.len().min(output.len()); - output[..copy_len].copy_from_slice(&decompressed[..copy_len]); - } else { - copy_chunk_to_output( - decompressed, - &mut output, - &chunk_offsets, - &chunk_dims, - &ds_dims, - &ds_strides, - &chunk_strides, - elem_size, - rank, - ); - } - } - - Ok(output) + None, + alloc_output, + |out| out.as_mut_slice(), + ) } /// Read a chunked dataset with caching support. /// /// On the first call, scans the chunk index (B-tree / fixed array / etc.) once /// and populates the cache's hash index. Subsequent calls skip the index scan -/// entirely. Decompressed chunk data is also cached with LRU eviction. +/// entirely. Decompressed chunk data is also cached with LRU eviction, when +/// the whole dataset fits in the cache. #[cfg(feature = "std")] #[allow(clippy::too_many_arguments)] pub fn read_chunked_data_cached( @@ -953,159 +1231,18 @@ pub fn read_chunked_data_cached( length_size: u8, cache: &ChunkCache, ) -> Result, FormatError> { - let (chunk_dimensions, version, addr_opt) = match layout { - DataLayout::Chunked { - chunk_dimensions, - version, - btree_address, - .. - } => (chunk_dimensions, *version, *btree_address), - _ => { - return Err(FormatError::ChunkedReadError( - "expected chunked layout".into(), - )); - } - }; - - let addr = addr_opt - .ok_or_else(|| FormatError::ChunkedReadError("no address for chunked layout".into()))?; - - check_chunk_element_size(layout, datatype, offset_size)?; - let elem_size = datatype.type_size() as usize; - let (rank, chunk_dims) = chunk_geometry(chunk_dimensions, version, dataspace, elem_size)?; - let ds_dims: Vec = dataspace.dimensions.iter().map(|&d| d as usize).collect(); - - // The per-file cache is shared across datasets (and threads); every - // lookup is keyed by this dataset's chunk-index address, so another - // dataset's index or chunks are never used for this read. - let chunks = cache.chunks_for(addr, rank, || { - list_chunks( - file_data, - layout, - dataspace, - elem_size, - offset_size, - length_size, - ) - .map(|(chunks, _)| chunks) - })?; - - // Assemble output - let total_bytes = checked_byte_len(dataspace.checked_num_elements()?, elem_size)?; - if total_bytes == 0 { - // Also keeps the stride products below in range: with a zero-sized - // dimension the total is 0 even if other dimensions are huge. - return Ok(Vec::new()); - } - let mut output = alloc_output(total_bytes)?; - - let mut ds_strides = vec![1usize; rank]; - for i in (0..rank.saturating_sub(1)).rev() { - ds_strides[i] = ds_strides[i + 1] * ds_dims[i + 1]; - } - - let mut chunk_strides = vec![1usize; rank]; - for i in (0..rank.saturating_sub(1)).rev() { - chunk_strides[i] = chunk_strides[i + 1] * chunk_dims[i + 1]; - } - - let chunk_total_bytes = checked_chunk_byte_len(&chunk_dims, elem_size)?; - - let mut place = |data: &[u8], chunk_info: &ChunkInfo| { - if rank == 0 { - let copy_len = data.len().min(output.len()); - output[..copy_len].copy_from_slice(&data[..copy_len]); - return; - } - let chunk_offsets: Vec = chunk_info - .offsets - .iter() - .take(rank) - .map(|&o| o as usize) - .collect(); - copy_chunk_to_output( - data, - &mut output, - &chunk_offsets, - &chunk_dims, - &ds_dims, - &ds_strides, - &chunk_strides, - elem_size, - rank, - ); - }; - let raw_bytes = |chunk_info: &ChunkInfo| -> Result<&[u8], FormatError> { - let c_addr = chunk_info.address as usize; - let size = chunk_info.chunk_size as usize; - ensure_len(file_data, c_addr, size)?; - Ok(&file_data[c_addr..c_addr + size]) - }; - - // Chunks stored as-is (no pipeline, or the filter mask says this chunk - // skipped every filter) are copied straight from the file bytes: they are - // already in memory, so routing them through a Vec and then an aligned - // cache buffer was two extra copies of the whole dataset for nothing. - let stored_raw = - |c: &ChunkInfo| pipeline.is_none_or(|pl| all_filters_skipped(pl, c.filter_mask)); - let mut misses: Vec<&ChunkInfo> = Vec::new(); - for chunk_info in &chunks { - if stored_raw(chunk_info) { - place(raw_bytes(chunk_info)?, chunk_info); - continue; - } - let coord: Vec = chunk_info.offsets.iter().take(rank).copied().collect(); - match cache.get_decompressed_in(addr, &coord) { - Some(cached) => place(&cached, chunk_info), - None => misses.push(chunk_info), - } - } - - // Decompress what the cache didn't have, a bounded batch at a time — in - // parallel with the `parallel` feature (this path, the one the facade - // uses, was sequential; only the uncached reader was parallel), unless the - // pool has one thread: then every reading thread would queue behind that - // one worker, so each decodes its own chunks instead. Chunks are - // cached only when the whole dataset fits: pushing a larger dataset - // through the cache just evicts each chunk moments after inserting it. - let cache_them = total_bytes <= cache.max_bytes(); - if let Some(pl) = pipeline { - let decode = |c: &&ChunkInfo| -> Result, FormatError> { - decompress_chunk_exact( - raw_bytes(c)?, - pl, - chunk_total_bytes, - elem_size as u32, - c.filter_mask, - &c.offsets, - ) - }; - for batch in misses.chunks(DECODE_BATCH) { - #[cfg(feature = "parallel")] - let decoded: Vec, FormatError>> = - if batch.len() >= 4 && parallel_read::pool_can_parallelise() { - use rayon::prelude::*; - batch.par_iter().map(decode).collect() - } else { - batch.iter().map(decode).collect() - }; - #[cfg(not(feature = "parallel"))] - let decoded: Vec, FormatError>> = batch.iter().map(decode).collect(); - - for (chunk_info, data) in batch.iter().zip(decoded) { - let data = data?; - if cache_them { - let coord: Vec = chunk_info.offsets.iter().take(rank).copied().collect(); - let cached = cache.put_decompressed_in(addr, coord, data); - place(&cached, chunk_info); - } else { - place(&data, chunk_info); - } - } - } - } - - Ok(output) + read_chunked_full( + file_data, + layout, + dataspace, + datatype, + pipeline, + offset_size, + length_size, + Some(cache), + alloc_output, + |out| out.as_mut_slice(), + ) } /// Sweep context passed into `read_chunked_data_sweep` to enable adaptive @@ -1486,6 +1623,7 @@ pub fn read_chunked_data_indexed( } /// Copy chunk data into the output buffer at the correct N-D position. +#[cfg(any(feature = "std", test))] #[allow(clippy::too_many_arguments)] fn copy_chunk_to_output( chunk_data: &[u8], @@ -1497,6 +1635,41 @@ fn copy_chunk_to_output( chunk_strides: &[usize], elem_size: usize, rank: usize, +) { + // SAFETY: `output` is exclusively borrowed, so no other write can + // overlap this one. + unsafe { + copy_chunk_into( + chunk_data, + &OutBuf::new(output), + chunk_offsets, + chunk_dims, + ds_dims, + ds_strides, + chunk_strides, + elem_size, + rank, + ) + } +} + +/// [`copy_chunk_to_output`] through an [`OutBuf`]: only the bytes of the +/// chunk's region (clipped to the dataset) are written. +/// +/// # Safety +/// +/// No other thread may be writing an overlapping region of `output`. +#[allow(clippy::too_many_arguments)] +unsafe fn copy_chunk_into( + chunk_data: &[u8], + output: &OutBuf<'_>, + chunk_offsets: &[usize], + chunk_dims: &[usize], + ds_dims: &[usize], + ds_strides: &[usize], + chunk_strides: &[usize], + elem_size: usize, + rank: usize, ) { // Row-copy approach: iterate over outer dimensions, memcpy the innermost // dimension in bulk. For 1-D data this is a single memcpy per chunk. @@ -1518,7 +1691,8 @@ fn copy_chunk_to_output( .is_some_and(|end| end <= output.len()) && src_bytes <= chunk_data.len() { - output[dst_start..dst_start + src_bytes].copy_from_slice(&chunk_data[..src_bytes]); + // SAFETY: in range (checked); exclusivity is the caller's contract. + unsafe { output.write(dst_start, &chunk_data[..src_bytes]) }; } return; } @@ -1620,8 +1794,8 @@ fn copy_chunk_to_output( .checked_add(row_bytes) .is_some_and(|end| end <= output.len()); if fits { - output[dst_start..dst_start + row_bytes] - .copy_from_slice(&chunk_data[src_start..src_start + row_bytes]); + // SAFETY: in range (checked); exclusivity is the caller's contract. + unsafe { output.write(dst_start, &chunk_data[src_start..src_start + row_bytes]) }; } } } @@ -2071,6 +2245,198 @@ mod tests { ); } + /// A 2-D `f32` dataset of `rows x cols` in `cr x cc` chunks, deflated + /// behind shuffle: (file bytes, chunk list, pipeline, expected output). + #[cfg(feature = "deflate")] + fn deflated_grid( + rows: usize, + cols: usize, + cr: usize, + cc: usize, + ) -> (Vec, Vec, FilterPipeline, Vec) { + use crate::filter_pipeline::{FILTER_DEFLATE, FILTER_SHUFFLE, FilterDescription}; + let pipeline = FilterPipeline { + version: 2, + filters: vec![ + FilterDescription { + filter_id: FILTER_SHUFFLE, + name: None, + flags: 0, + client_data: vec![4], + }, + FilterDescription { + filter_id: FILTER_DEFLATE, + name: None, + flags: 0, + client_data: vec![1], + }, + ], + }; + let value = |r: usize, c: usize| (r * cols + c) as f32 * 0.5; + let expected: Vec = (0..rows * cols) + .flat_map(|i| value(i / cols, i % cols).to_le_bytes()) + .collect(); + let (mut file, mut chunks) = (vec![0u8; 8], Vec::new()); + for r0 in (0..rows).step_by(cr) { + for c0 in (0..cols).step_by(cc) { + // Full-size chunks; the edge padding is zeros. + let chunk: Vec = (0..cr * cc) + .flat_map(|i| { + let (r, c) = (r0 + i / cc, c0 + i % cc); + let v = if r < rows && c < cols { + value(r, c) + } else { + 0.0 + }; + v.to_le_bytes() + }) + .collect(); + let stored = crate::filters::compress_chunk(&chunk, &pipeline, 4).unwrap(); + chunks.push(ChunkInfo { + chunk_size: stored.len() as u32, + filter_mask: 0, + offsets: vec![r0 as u64, c0 as u64, 0], + address: file.len() as u64, + }); + file.extend(stored); + } + } + (file, chunks, pipeline, expected) + } + + #[cfg(feature = "deflate")] + fn fill( + file: &[u8], + chunks: &[ChunkInfo], + pipeline: &FilterPipeline, + dims: [usize; 2], + chunk: [usize; 2], + ) -> Result, FormatError> { + let placer = ChunkPlacer::new(&chunk, &dims, 4); + let mut out = vec![0u8; dims[0] * dims[1] * 4]; + fill_from_chunks( + file, + chunks, + Some(pipeline), + &placer, + chunk[0] * chunk[1] * 4, + None, + &mut out, + )?; + Ok(out) + } + + /// Chunks decoded straight into the output (in parallel with the + /// `parallel` feature, partial edge chunks included) land where the + /// row-by-row reference puts them. + #[test] + #[cfg(feature = "deflate")] + fn chunks_fill_the_output_in_place() { + for (dims, chunk) in [ + ([100, 70], [16, 32]), + ([64, 64], [16, 16]), + ([5, 3], [8, 8]), + ] { + let (file, chunks, pipeline, expected) = + deflated_grid(dims[0], dims[1], chunk[0], chunk[1]); + for _ in 0..4 { + assert_eq!( + fill(&file, &chunks, &pipeline, dims, chunk).unwrap(), + expected + ); + } + // Any chunk order gives the same output. + let mut reversed = chunks.clone(); + reversed.reverse(); + assert_eq!( + fill(&file, &reversed, &pipeline, dims, chunk).unwrap(), + expected + ); + } + } + + /// Of several corrupt chunks, the error names the first in chunk order, + /// however the chunks were shared out between threads. + #[test] + #[cfg(feature = "deflate")] + fn first_corrupt_chunk_is_the_error() { + let (dims, chunk) = ([128, 64], [8, 16]); + let (mut file, mut chunks, pipeline, _) = + deflated_grid(dims[0], dims[1], chunk[0], chunk[1]); + // Valid streams that decode short: an error naming the chunk. + let short = crate::filters::compress_chunk(&[1u8; 64], &pipeline, 4).unwrap(); + for bad in [5usize, 11, 40] { + chunks[bad].address = file.len() as u64; + chunks[bad].chunk_size = short.len() as u32; + file.extend_from_slice(&short); + } + for _ in 0..20 { + let err = fill(&file, &chunks, &pipeline, dims, chunk).unwrap_err(); + let want = format!("{:?}", chunks[5].offsets); + assert!(err.to_string().contains(&want), "{err} (want {want})"); + } + } + + /// Only chunks on the grid, each at a distinct place, may be placed + /// concurrently; anything else is read one chunk at a time. + #[test] + #[cfg(feature = "parallel")] + fn regions_disjoint_only_for_distinct_grid_chunks() { + let placer = ChunkPlacer::new(&[10, 10], &[25, 30], 4); + let chunk = |r: u64, c: u64| ChunkInfo { + chunk_size: 0, + filter_mask: 0, + offsets: vec![r, c, 0], + address: 0, + }; + let grid: Vec = (0..3) + .flat_map(|r| (0..3).map(move |c| chunk(r * 10, c * 10))) + .collect(); + assert!(placer.regions_disjoint(&grid)); + // Chunks wholly outside the dataset write nothing. + let mut outside = grid.clone(); + outside.push(chunk(30, 0)); + outside.push(chunk(30, 0)); + assert!(placer.regions_disjoint(&outside)); + let mut duplicate = grid.clone(); + duplicate.push(chunk(10, 20)); + assert!(!placer.regions_disjoint(&duplicate)); + let mut off_grid = grid.clone(); + off_grid[4] = chunk(15, 10); + assert!(!placer.regions_disjoint(&off_grid)); + let mut short = grid; + short[0].offsets.truncate(1); + assert!(!placer.regions_disjoint(&short)); + } + + /// A duplicated chunk in a corrupt index is not placed from two threads + /// at once: the read goes one chunk at a time, as before. + #[test] + #[cfg(feature = "deflate")] + fn duplicate_chunks_are_read_in_order() { + let (dims, chunk) = ([64, 64], [16, 16]); + let (file, mut chunks, pipeline, expected) = + deflated_grid(dims[0], dims[1], chunk[0], chunk[1]); + // A second entry at chunk 3's place, holding chunk 7's data: the + // later entry wins, as a sequential read has it. + let mut dup = chunks[7].clone(); + dup.offsets = chunks[3].offsets.clone(); + chunks.push(dup); + let got = fill(&file, &chunks, &pipeline, dims, chunk).unwrap(); + let chunk7_in_3: Vec = { + let mut e = expected.clone(); + let (r3, c3) = (chunks[3].offsets[0] as usize, chunks[3].offsets[1] as usize); + let (r7, c7) = (chunks[7].offsets[0] as usize, chunks[7].offsets[1] as usize); + for r in 0..16 { + let src = ((r7 + r) * 64 + c7) * 4; + let dst = ((r3 + r) * 64 + c3) * 4; + e[dst..dst + 64].copy_from_slice(&expected[src..src + 64]); + } + e + }; + assert_eq!(got, chunk7_in_3); + } + #[test] fn copy_chunk_to_output_1d_rejects_overflowing_offset_without_panicking() { // Found by fuzzing: `global_start * elem_size` overflowed for a @@ -2481,6 +2847,68 @@ mod tests { use crate::chunk_cache::ChunkCache; + /// A chunk stored unfiltered in a filtered dataset (every filter-mask + /// bit set) must still hold the whole chunk. The cached reader placed a + /// short one and read its missing rows as zeros; both readers refuse it. + #[test] + fn short_unfiltered_chunk_of_a_filtered_dataset_is_an_error() { + use crate::filter_pipeline::{FILTER_SHUFFLE, FilterDescription}; + let mut file_data = vec![0u8; 0x2000]; + let mut infos = Vec::new(); + for k in 0..3u64 { + let address = 0x1000 + k as usize * 80; + for i in 0..10u64 { + let at = address + i as usize * 8; + file_data[at..at + 8].copy_from_slice(&((k * 10 + i) as f64).to_le_bytes()); + } + infos.push(ChunkInfo { + chunk_size: if k == 1 { 40 } else { 80 }, + filter_mask: 1, + offsets: vec![k * 10, 0], + address: address as u64, + }); + } + let btree = build_chunk_btree_leaf(&infos, 2, 8); + file_data[0x100..0x100 + btree.len()].copy_from_slice(&btree); + let layout = DataLayout::Chunked { + chunk_dimensions: vec![10, 8], + btree_address: Some(0x100), + version: 3, + chunk_index_type: None, + single_chunk_filtered_size: None, + single_chunk_filter_mask: None, + dont_filter_partial_edge_chunks: false, + }; + let dataspace = simple_space(vec![30]); + let pipeline = FilterPipeline { + version: 2, + filters: vec![FilterDescription { + filter_id: FILTER_SHUFFLE, + name: None, + flags: 0, + client_data: vec![8], + }], + }; + let dt = make_f64_type(); + let cache = ChunkCache::new(); + for result in [ + read_chunked_data(&file_data, &layout, &dataspace, &dt, Some(&pipeline), 8, 8), + read_chunked_data_cached( + &file_data, + &layout, + &dataspace, + &dt, + Some(&pipeline), + 8, + 8, + &cache, + ), + ] { + let err = result.unwrap_err().to_string(); + assert!(err.contains("[10, 0]") && err.contains("40"), "{err}"); + } + } + #[test] fn cached_read_populates_index_and_returns_correct_data() { let values: Vec = (0..20).map(|i| i as f64).collect(); diff --git a/crates/clawhdf5-format/src/parallel_read.rs b/crates/clawhdf5-format/src/parallel_read.rs index 1b940c6..61d06cd 100644 --- a/crates/clawhdf5-format/src/parallel_read.rs +++ b/crates/clawhdf5-format/src/parallel_read.rs @@ -41,6 +41,126 @@ pub fn pool_can_parallelise() -> bool { rayon::current_num_threads() > 1 } +/// How many rayon workers [`run_with_helpers`] should ask to help with +/// `items` work items, given that the calling thread works too: the pool's +/// other threads (all of them when the caller is not one), at most one per +/// item beyond the caller's first. +pub(crate) fn helper_count(items: usize) -> usize { + let pool = rayon::current_num_threads(); + let others = if rayon::current_thread_index().is_some() { + pool.saturating_sub(1) + } else { + pool + }; + others.min(items.saturating_sub(1)) +} + +/// Run `body` on the calling thread and on up to `helpers` rayon workers at +/// once, returning when the caller's call has finished and every worker that +/// started one has too. `body` shares its work out itself (typically by +/// claiming items from an atomic counter until none are left). +/// +/// The caller never waits for a worker to *become* free: helpers are queued +/// on the pool, and one that only gets to run after the caller has finished +/// returns without calling `body`. So a busy or small pool can only fail to +/// speed a read up, never hold it back — with `par_iter`, the calling thread +/// (not a pool worker) handed all the work to the pool and slept, and N +/// threads reading through a 2-worker pool decoded on 2 cores. +/// +/// A panic in `body`, on any thread, is resumed on the caller once every +/// helper that started has stopped. +pub(crate) fn run_with_helpers(helpers: usize, body: &(dyn Fn() + Sync)) { + use std::panic::{AssertUnwindSafe, catch_unwind, resume_unwind}; + use std::sync::{Arc, Condvar, Mutex, PoisonError}; + + if helpers == 0 { + body(); + return; + } + + type Body = dyn Fn() + Sync + 'static; + struct Shared { + /// `body`, its lifetime erased. Only dereferenced by a helper that + /// registered in `state` while it was open (see below). + body: *const Body, + /// (closed, helpers inside `body`). + state: Mutex<(bool, usize)>, + idle: Condvar, + panic: Mutex>>, + } + // SAFETY: `body` points to a `Sync` closure, so calling it from other + // threads is allowed; the pointer is only used under the protocol below, + // which keeps it from outliving the closure. + unsafe impl Send for Shared {} + unsafe impl Sync for Shared {} + + fn help(shared: &Shared) { + { + let mut state = shared.state.lock().unwrap_or_else(PoisonError::into_inner); + if state.0 { + return; + } + state.1 += 1; + } + // SAFETY: registered while open, so the caller of `run_with_helpers` + // is still inside it (it closes, then waits until no helper is + // registered, before returning), and `body` is alive. + let body = unsafe { &*shared.body }; + if let Err(payload) = catch_unwind(AssertUnwindSafe(body)) { + shared + .panic + .lock() + .unwrap_or_else(PoisonError::into_inner) + .get_or_insert(payload); + } + let mut state = shared.state.lock().unwrap_or_else(PoisonError::into_inner); + state.1 -= 1; + if state.1 == 0 { + shared.idle.notify_all(); + } + } + + let body_ptr: *const (dyn Fn() + Sync + '_) = body; + // SAFETY: only the lifetime changes (same fat-pointer layout). The + // pointer is dereferenced only while this function is running: see + // `help` and the wait below. + let body_ptr: *const Body = unsafe { core::mem::transmute(body_ptr) }; + let shared = Arc::new(Shared { + body: body_ptr, + state: Mutex::new((false, 0)), + idle: Condvar::new(), + panic: Mutex::new(None), + }); + for _ in 0..helpers { + let shared = Arc::clone(&shared); + rayon::spawn(move || help(&shared)); + } + let caller = catch_unwind(AssertUnwindSafe(body)); + { + // Close, then wait for the helpers inside `body`; later ones return + // at once. This must happen even if `body` panicked on this thread. + let mut state = shared.state.lock().unwrap_or_else(PoisonError::into_inner); + state.0 = true; + while state.1 > 0 { + state = shared + .idle + .wait(state) + .unwrap_or_else(PoisonError::into_inner); + } + } + if let Err(payload) = caller { + resume_unwind(payload); + } + let helper_panic = shared + .panic + .lock() + .unwrap_or_else(PoisonError::into_inner) + .take(); + if let Some(payload) = helper_panic { + resume_unwind(payload); + } +} + /// Decompress chunks in parallel using lane-partitioned assignment. /// /// Instead of naive `par_iter`, chunks are deterministically assigned to lanes @@ -258,6 +378,60 @@ mod tests { (file, infos) } + /// Every item is processed exactly once, whatever mix of caller and + /// helpers ends up doing it. + #[test] + fn run_with_helpers_shares_all_work() { + use core::sync::atomic::{AtomicUsize, Ordering}; + for helpers in [0, 1, 3, 16] { + let n = 1000; + let next = AtomicUsize::new(0); + let done: Vec = (0..n).map(|_| AtomicUsize::new(0)).collect(); + run_with_helpers(helpers, &|| { + loop { + let i = next.fetch_add(1, Ordering::Relaxed); + if i >= n { + break; + } + done[i].fetch_add(1, Ordering::Relaxed); + } + }); + assert!(done.iter().all(|d| d.load(Ordering::Relaxed) == 1)); + } + } + + /// A panic in the shared body reaches the caller whichever thread it + /// happened on, and only after the helpers inside the body have left it + /// (they borrow the caller's stack). + #[test] + fn run_with_helpers_propagates_panics() { + use core::sync::atomic::{AtomicUsize, Ordering}; + use std::panic::{AssertUnwindSafe, catch_unwind}; + let caller = std::thread::current().id(); + for panic_on_caller in [true, false] { + let inside = AtomicUsize::new(0); + let calls = AtomicUsize::new(0); + let result = catch_unwind(AssertUnwindSafe(|| { + run_with_helpers(4, &|| { + inside.fetch_add(1, Ordering::SeqCst); + calls.fetch_add(1, Ordering::SeqCst); + let on_caller = std::thread::current().id() == caller; + std::thread::sleep(std::time::Duration::from_millis(20)); + inside.fetch_sub(1, Ordering::SeqCst); + if on_caller == panic_on_caller { + panic!("boom"); + } + }); + })); + // A helper may never have run (the pool was slow to start it), + // in which case nothing panicked when `panic_on_caller` is false. + if panic_on_caller || calls.load(Ordering::SeqCst) > 1 { + assert!(result.is_err()); + } + assert_eq!(inside.load(Ordering::SeqCst), 0); + } + } + /// Every parallel decoder refuses a chunk that decodes short, naming it. #[test] fn short_decoded_chunk_is_an_error() { diff --git a/crates/clawhdf5/tests/busy_decode_pool.rs b/crates/clawhdf5/tests/busy_decode_pool.rs new file mode 100644 index 0000000..d53eeab --- /dev/null +++ b/crates/clawhdf5/tests/busy_decode_pool.rs @@ -0,0 +1,112 @@ +//! Full reads of chunked datasets must not wait for a busy rayon pool of +//! any size. +//! +//! A full read handed its chunks to the rayon pool (`par_iter`) and the +//! calling thread — not a pool worker — slept until the pool had decoded +//! them. With a small pool (2-4 threads) and more reading threads than +//! workers, every reader queued behind the same few workers +//! (`docs/known-issues.md`, "Concurrent and contiguous read performance"). +//! Now the calling thread decodes too, and pool workers only help when they +//! are free. The test keeps both workers of a two-thread pool busy and +//! requires reads to finish anyway, with the right values. +//! +//! One test in its own binary: it configures the process-wide rayon pool. + +#![cfg(feature = "parallel")] + +use std::sync::mpsc; +use std::time::Duration; + +use clawhdf5::{File, FileBuilder}; + +const N: usize = 4096; // 64 chunks of 64 elements + +fn values() -> Vec { + (0..N).map(|i| i as f64 * 0.25 - 7.0).collect() +} + +fn build() -> File { + let mut b = FileBuilder::new(); + b.create_dataset("data") + .with_f64_data(&values()) + .with_shape(&[N as u64]) + .with_chunks(&[64]) + .with_deflate(1) + .with_provenance("test-suite", "2026-09-26T00:00:00Z", None); + File::from_bytes(b.finish().unwrap()).unwrap() +} + +/// Run `f` on a fresh thread; `None` if it has not finished within `limit`. +fn finishes_within( + limit: Duration, + f: impl FnOnce() -> T + Send + 'static, +) -> Option { + let (tx, rx) = mpsc::channel(); + std::thread::spawn(move || { + let _ = tx.send(f()); + }); + rx.recv_timeout(limit).ok() +} + +#[test] +fn full_reads_do_not_wait_for_a_busy_small_pool() { + const WORKERS: usize = 2; + rayon::ThreadPoolBuilder::new() + .num_threads(WORKERS) + .build_global() + .expect("this test binary configures the global pool first"); + + // Built first: the writer compresses on the pool too. + let file = std::sync::Arc::new(build()); + + // Occupy every worker of the pool until the reads are done. + let (started_tx, started_rx) = mpsc::channel(); + let (release_tx, release_rx) = mpsc::channel::<()>(); + let release_rx = std::sync::Arc::new(std::sync::Mutex::new(release_rx)); + for _ in 0..WORKERS { + let (started_tx, release_rx) = (started_tx.clone(), release_rx.clone()); + rayon::spawn(move || { + started_tx.send(()).unwrap(); + let _ = release_rx.lock().unwrap().recv(); + }); + } + for _ in 0..WORKERS { + started_rx.recv().unwrap(); + } + + let limit = Duration::from_secs(20); + // Several readers at once, as in the `concurrent_read` benchmark: the + // cached full read (`read_*`), the typed one and the uncached reader + // behind `verify_provenance`. + let readers: Vec<_> = (0..4) + .map(|_| { + let file = std::sync::Arc::clone(&file); + std::thread::spawn(move || { + finishes_within(limit, move || { + let ds = file.dataset("data").unwrap(); + ( + ds.read_f64().unwrap(), + ds.read_f32().unwrap(), + ds.verify_provenance().unwrap(), + ) + }) + }) + }) + .collect(); + let results: Vec<_> = readers.into_iter().map(|h| h.join().unwrap()).collect(); + // Free the workers before asserting, so a failure does not hang the + // blocked reader threads forever. + for _ in 0..WORKERS { + release_tx.send(()).unwrap(); + } + + let want = values(); + let want_f32: Vec = want.iter().map(|&v| v as f32).collect(); + for result in results { + let (f64s, f32s, verified) = + result.expect("a full read waited for the busy two-thread rayon pool"); + assert_eq!(f64s, want); + assert_eq!(f32s, want_f32); + assert_eq!(verified, clawhdf5::provenance::VerifyResult::Ok); + } +}