Packed fingerprint corpora¶
smsd.corpus.FingerprintCorpus stores binary fingerprints as fixed-width packed
rows in a read-only, memory-mapped file. NumPy is required for this API and for
smsd.clustering; importing the modules does not import NumPy.
For a task-based walk-through of building, reopening and searching an index, see the packed fingerprint index chapter of the user guide.
from pathlib import Path
from tempfile import TemporaryDirectory
from smsd.corpus import FingerprintCorpus
from smsd.clustering import butina, maxmin
with TemporaryDirectory() as directory:
path = Path(directory) / "library.fp"
with FingerprintCorpus.create(path, [[1], [1], [2]], fp_size=8) as corpus:
assert len(corpus) == 3
assert corpus.fingerprint(0) == [1]
assert list(corpus.similarities([1], threshold=1.0)) == [(0, 1.0), (1, 1.0)]
assert list(corpus.screen_many([[1], [2]], threshold=1.0)) == [
(0, 0, 1.0), (0, 1, 1.0), (1, 2, 1.0)
]
assert list(corpus.threshold_pairs(1.0)) == [(0, 1, 1.0)]
assert list(corpus.top_k_pairs(1)) == [
(0, 1, 1.0), (1, 0, 1.0), (2, 0, 0.0)
]
assert butina(corpus, threshold=1.0) == [(0, 1), (2,)]
assert maxmin(corpus, count=3, seed=0) == [0, 2, 1]
chunks = list(corpus.iter_packed_chunks(chunk_size=2))
assert chunks[0][0] == 0
assert chunks[0][1].tolist() == [[2], [2]]
assert corpus.cuda_calls == 0 # The default backend is CPU.
# Chunk arrays own their data and remain usable after the mapping is closed.
assert chunks[1][1].tolist() == [[4]]
with FingerprintCorpus(path) as reopened:
assert reopened.fingerprint(2) == [2]
from_fingerprints is an alias for create. The fingerprint input can be a
generator. Construction publishes a complete file atomically, preserves an
existing destination by default and removes its temporary file after an input
failure. Use overwrite=True for replacement. The containing filesystem must
support hard links for the default create mode.
To create packed rows directly from molecule inputs, use from_smiles. Queries
must use the same fingerprint specification:
from pathlib import Path
from tempfile import TemporaryDirectory
import smsd
from smsd.corpus import FingerprintCorpus
molecules = ["CCO", "CCN", "c1ccccc1", "CC(=O)O"]
options = smsd.FingerprintOptions.ecfp(radius=2, size=2048)
with TemporaryDirectory() as directory:
with FingerprintCorpus.from_smiles(
Path(directory) / "molecules.fp", iter(molecules), options=options,
chunk_size=1024
) as corpus:
query = smsd.fingerprint("CCO", options=options)
for index, score in corpus.similarities(query, threshold=0.3):
print(index, score)
The file stores its bit width and row count, rather than the chemistry
specification. Keep that specification with the application. Count fingerprints
are rejected; duplicate binary positions are folded. Two empty fingerprints
score 1, and an empty fingerprint compared with a non-empty one scores 0, as
with smsd.tanimoto_coefficient. Thresholds are inclusive and must be finite
values from 0 to 1.
CPU and CUDA backends¶
screen_many, similarities, threshold_pairs, top_k_pairs, butina and
maxmin accept the same explicit backend argument:
| Backend | Behaviour |
|---|---|
"cpu" |
Default. NumPy computes packed intersections; no CUDA provider is requested. |
"cuda" |
Requires a compatible smsd-cuda package and device. Unavailable providers or unsupported widths raise. |
"auto" |
Uses the CUDA provider when available and supported; otherwise uses CPU comparisons. |
The core and CUDA package versions must match. Packed CUDA supports fingerprint
widths up to 65,536 bits and a little-endian host. Wider packed fingerprints
remain supported on the CPU; auto uses that CPU path. Errors during a device
call, result validation or resident cleanup propagate. They do not silently
retry the operation on the CPU. These packed operations do not use the native
TargetCorpus operation crossovers.
A CUDA-capable installation can request the device explicitly:
from pathlib import Path
from tempfile import TemporaryDirectory
from smsd.corpus import FingerprintCorpus
from smsd.clustering import butina, maxmin
with TemporaryDirectory() as directory:
with FingerprintCorpus.create(
Path(directory) / "bits.fp", [[1], [1, 2], [3]], fp_size=8
) as corpus:
before = corpus.cuda_calls
hits = list(corpus.screen_many([[1], [3]], threshold=0.5, backend="cuda"))
assert hits == [(0, 0, 1.0), (0, 1, 0.5), (1, 2, 1.0)]
assert corpus.cuda_calls > before
pairs = list(corpus.threshold_pairs(0.5, backend="cuda"))
neighbours = list(corpus.top_k_pairs(1, backend="cuda"))
clusters = butina(corpus, threshold=0.5, backend="cuda")
selected = maxmin(corpus, count=2, backend="cuda")
CUDA computes exact integer intersections. screen_many applies integer cutoff
tables on device, chosen to reproduce the CPU's binary64 threshold comparison;
the host forms the returned Tanimoto scores from integer counts. top_k_pairs
performs threshold selection and per-tile top-k on device, then merges those
lists on the CPU with the same descending score and ascending target-index tie
break. Pair searches and clustering use the same exact scores. GPU acceleration
does not make those algorithms approximate.
Butina centre selection, count updates and assignment bookkeeping remain on the
CPU. MaxMin's minimum-distance updates and next-row choice also remain on the
CPU. Both have deterministic lower-index tie breaks. The cuda_calls property
is a cumulative count of completed packed CUDA calls, including clustering
passes. A positive increase is evidence that this adapter completed CUDA work;
it is not a count of queries, hits or GPU kernels, and does not establish a
speedup or the fraction of an application executed on the GPU.
Ordering, work and memory¶
screen_many accepts at most 32 packed queries and consumes at most one extra
query input to detect an excess before reading targets. It yields target-chunk
order, then query index and target index. This is not a global query-major or
similarity ranking. similarities yields target input order for one query.
CUDA tiles contain at most 1,048,576 query/target pairs. A conservative 64 MiB
host-work estimate also limits each tile, so the effective target chunk can be
smaller than the requested chunk_size. CPU packed comparison temporaries are
similarly bounded to approximately 64 MiB; chunk_size can lower the bound.
These are array/work estimates, not hard caps on process memory. Mapped pages,
query storage, Python objects and application results are additional costs.
Molecule construction and parsing memory are also input-dependent.
The native resident workspace separately retains at most 256 MiB of device
scratch. That cap is not a total GPU-memory limit: resident target fingerprints,
other temporary allocations, CUDA runtime resources and host pinned staging are
separate. The packed adapter creates and releases a target resident handle for
each tile; an open memory-mapped file does not imply that the entire corpus is
resident on the GPU or that targets persist on device between calls. The native
TargetCorpus API has its own retained-library lifetime.
top_k_pairs yields directed neighbours, ordered by descending similarity and
then target index, with at most k heap entries in the global merge.
threshold_pairs yields each undirected pair once in index order. Both exclude
self pairs by default. Neither creates a dense all-against-all matrix. Retaining
every yielded threshold pair can nevertheless require quadratic application
memory, and retaining all top-k results can require O(n * k) memory.
The exact all-against-all scans take quadratic comparison work. Butina keeps
linear auxiliary arrays and returns tuples with the centre first. Its default
uses initial neighbour counts; reordering=True updates those counts after
each cluster. A member must meet the similarity threshold with its centre,
but need not meet it with every other member. MaxMin keeps a linear array of
minimum Tanimoto distances and uses a deterministic starting row; comparison
work is O(n * selected_count).
The format has a 24-byte little-endian header and ceil(fp_size / 8) bytes per
row. The supported width is 1 to 16,777,216 bits. A 2048-bit corpus requires
256 bytes per row plus the header. The format and bounded tiles do not establish
validation at a particular corpus size or any throughput result.
Keep an open corpus file unchanged, and close it with a context manager.
Concurrent reads copy their mapped rows under a shared lifetime lock. close()
waits for an active copy before unmapping; repeated closes are harmless.
Computation on copied rows and already yielded results remains valid. A running
iterator raises ValueError when it next requests a mapped row after close.
Closing does not cancel a whole search. Explicit lifecycle calls are available
when a context manager is unsuitable:
from pathlib import Path
from tempfile import TemporaryDirectory
from smsd.corpus import FingerprintCorpus
with TemporaryDirectory() as directory:
corpus = FingerprintCorpus.create(Path(directory) / "bits.fp", [[1]], fp_size=8)
try:
assert corpus.fingerprint(0) == [1]
finally:
corpus.close()
corpus.close()
Streamed native screening¶
stream_screen_smiles screens a target iterable with the public native corpus
API and releases each target chunk before loading the next. It uses the fixed
ECFP radius 2, 2048-bit specification of TargetCorpus.screen. It requires no
NumPy installation. This native graph path is separate from the packed backend
choices above.
from smsd.corpus import stream_screen_smiles
targets = iter(["CCO", "CCN", "c1ccccc1", "CCO"])
for index in stream_screen_smiles("CCO", targets, 0.5, chunk_size=4096):
print(index)
Each chunk holds at most chunk_size target inputs, their prepared native
graphs and one chunk's hit indices; the maximum accepted chunk is 1,048,576
targets. Graph memory depends on the input molecules. Hit indices retain their
global input offsets. An invalid later molecule raises after any earlier
chunks' hits have been yielded.
The native screen uses the optional smsd-cuda backend when eligible and
retains its ordinary CPU fallback. The operation's crossover applies to each
chunk:
a small chunk may run on the CPU even when a GPU is available. This wrapper
does not force GPU execution. For repeated queries over one retained library,
use TargetCorpus directly so its resident device allocation can be reused.
For several query inputs, stream_screen_many_smiles uses bounded query and
target tiles through the same native operation:
from smsd.corpus import stream_screen_many_smiles
molecules = ["CCO", "CCN", "c1ccccc1", "CC(=O)O"]
for query_index, target_index in stream_screen_many_smiles(
iter(["CCO", "CCN"]), iter(molecules), 0.5,
query_chunk_size=64, target_chunk_size=4096, max_queries=4096
):
print(query_index, target_index)
The query inventory is limited by max_queries; one extra query input may be
consumed to reject an excess before reading targets. Queries are prepared once,
and the target iterable is consumed once. Native result storage covers at most
one query tile by one target tile. All three controls must be positive and at
most 1,048,576, and the query/target tile product cannot exceed 1,048,576.
These are bounds on molecule counts and result pairs; graph memory depends on
the input molecules. Hits follow target-tile order, then query index, then
target index. They are not globally query-major or ranked by similarity.
The packed corpus's 32-query limit does not apply to this native tiled wrapper.
Empty queries yield no hits without consuming targets. Invalid later targets
can raise after earlier tiles have yielded hits.
smsd.warm_start(corpus, query, threshold=0.5) discovers available accelerators,
prewarms that real corpus if needed and returns its ordinary screen hit indices.
The corpus retains any prepared fingerprints and resident device allocation.
The same dispatch rules apply; a successful call alone does not establish that
a GPU kernel ran.