Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
965 changes: 305 additions & 660 deletions benchmark/results.json

Large diffs are not rendered by default.

167 changes: 0 additions & 167 deletions benchmark/suite/layout.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,105 +2,20 @@

from typing import Any

import numpy as np
from zdisamar import rtm

from . import cases, config

INLINE_INTEGRATION_SAMPLE_COUNT = 5
WAVELENGTH_SAMPLING_ROW_BYTES = 200
OLD_WAVELENGTH_SAMPLING_ROW_BYTES = 65_592
OWNED_WAVELENGTH_SAMPLING_HEADER_BYTES = 48
INTEGRATION_KERNEL_BYTES = 32_784
SIDE_SAMPLE_BYTES = 16
MIN_PARALLEL_WAVELENGTH_SAMPLE_COUNT = 64
MAX_WORKERS = 64
COLLISION_COMPLEX_PROFILE_CACHE_OLD_BYTES = 2_072
COLLISION_COMPLEX_PROFILE_CACHE_BYTES = 1_048
COLLISION_COMPLEX_PROFILE_CACHE_CAPACITY = 64
COLLISION_COMPLEX_PROFILE_NODE_COUNT = 47
SPLINE_ENDPOINT_SECANT_SAMPLE_SCRATCH_BYTES = 10_240

INTEGRATION_MODE_LABELS = {
0: "auto",
1: "explicit_hr_grid",
2: "disamar_hr_grid",
3: "adaptive",
}
CHANNEL_LABELS = {
0: "radiance",
1: "irradiance",
}


def memory_layout_diagnostics() -> dict[str, Any]:

return {
"instrument_sampling": instrument_sampling_layout(),
"optical_accumulation": optical_accumulation_layout(),
}
Comment thread
bout3fiddy marked this conversation as resolved.


def instrument_sampling_layout() -> dict[str, Any]:

case = cases.forward_case()
nominal_wavelengths_nm = rtm.nominal_wavelengths(case)
response = rtm.instrument_response(case, nominal_wavelengths_nm)
rows = response.to_rows()
kernel_rows = [row for row in rows if int(row["sample_index"]) == 0]
support_counts = np.asarray([row["support_count"] for row in kernel_rows], dtype=np.int64)
nominal_sample_count = len(nominal_wavelengths_nm)
kernel_count = int(support_counts.size)
side_sample_count = int(support_counts[support_counts > INLINE_INTEGRATION_SAMPLE_COUNT].sum())
side_kernel_count = int(np.count_nonzero(support_counts > INLINE_INTEGRATION_SAMPLE_COUNT))
current_row_payload_bytes = nominal_sample_count * WAVELENGTH_SAMPLING_ROW_BYTES
current_side_storage_bytes = side_sample_count * SIDE_SAMPLE_BYTES
current_owned_plan_bytes = (
OWNED_WAVELENGTH_SAMPLING_HEADER_BYTES
+ current_row_payload_bytes
+ current_side_storage_bytes
)
previous_row_payload_bytes = nominal_sample_count * OLD_WAVELENGTH_SAMPLING_ROW_BYTES
saved_row_payload_bytes = previous_row_payload_bytes - (
current_row_payload_bytes + current_side_storage_bytes
)
sampling_worker_count = preferred_sampling_worker_count(nominal_sample_count)

return {
"boundary": "post-timing diagnostic; excluded from benchmark wall time and peak RSS",
"source": "instrument_response_sampling over benchmark forward_case nominal grid",
"nominal_sample_count": nominal_sample_count,
"kernel_count": kernel_count,
"support_sample_count": len(rows),
"support_count_stats": int_stats(support_counts),
"support_count_histogram": int_histogram(support_counts),
"integration_mode_histogram": labeled_histogram(
np.asarray([row["integration_mode"] for row in kernel_rows], dtype=np.int64),
INTEGRATION_MODE_LABELS,
),
"channel_stats": channel_stats(kernel_rows),
"inline_threshold_samples": INLINE_INTEGRATION_SAMPLE_COUNT,
"inline_kernel_count": int(kernel_count - side_kernel_count),
"side_kernel_count": side_kernel_count,
"side_sample_count": side_sample_count,
"estimated_sampling_worker_count": sampling_worker_count,
"integration_kernel_scratch_bytes_per_worker": INTEGRATION_KERNEL_BYTES,
"integration_kernel_scratch_bytes_at_worker_count": (
sampling_worker_count * INTEGRATION_KERNEL_BYTES
),
"current_owned_wavelength_sampling_bytes": current_owned_plan_bytes,
"current_row_payload_bytes": current_row_payload_bytes,
"current_side_storage_bytes": current_side_storage_bytes,
"previous_full_kernel_row_payload_bytes": previous_row_payload_bytes,
"estimated_row_payload_saved_bytes": saved_row_payload_bytes,
"estimated_row_payload_saved_mib": mib(saved_row_payload_bytes),
"estimated_row_payload_reduction_ratio": ratio(
saved_row_payload_bytes,
previous_row_payload_bytes,
),
}


def optical_accumulation_layout() -> dict[str, Any]:

cache_saved_bytes = (
Expand Down Expand Up @@ -131,88 +46,6 @@ def optical_accumulation_layout() -> dict[str, Any]:
}


def channel_stats(kernel_rows: list[dict[str, Any]]) -> dict[str, Any]:

result = {}

for channel_code, label in CHANNEL_LABELS.items():
rows = [row for row in kernel_rows if int(row["channel"]) == channel_code]
counts = np.asarray([row["support_count"] for row in rows], dtype=np.int64)
side_counts = counts[counts > INLINE_INTEGRATION_SAMPLE_COUNT]
result[label] = {
"kernel_count": int(counts.size),
"support_count_stats": int_stats(counts),
"support_count_histogram": int_histogram(counts),
"side_kernel_count": int(side_counts.size),
"side_sample_count": int(side_counts.sum()),
}

return result


def preferred_sampling_worker_count(sample_count: int) -> int:

if sample_count < MIN_PARALLEL_WAVELENGTH_SAMPLE_COUNT:
return 1

count_from_work = max(1, sample_count // MIN_PARALLEL_WAVELENGTH_SAMPLE_COUNT)
worker_cap = config.EFFECTIVE_WORKER_CAP or 1

return min(MAX_WORKERS, worker_cap, count_from_work)


def int_stats(values: np.ndarray) -> dict[str, float | int]:

if values.size == 0:
return {
"count": 0,
"min": 0,
"median": 0.0,
"mean": 0.0,
"p90": 0.0,
"max": 0,
}

return {
"count": int(values.size),
"min": int(np.min(values)),
"median": float(np.median(values)),
"mean": float(np.mean(values)),
"p90": float(np.percentile(values, 90.0)),
"max": int(np.max(values)),
}


def int_histogram(values: np.ndarray) -> dict[str, int]:

if values.size == 0:
return {}

unique, counts = np.unique(values.astype(np.int64), return_counts=True)

return {str(int(value)): int(count) for value, count in zip(unique, counts, strict=True)}


def labeled_histogram(values: np.ndarray, labels: dict[int, str]) -> dict[str, int]:

histogram = int_histogram(values)

return {labels.get(int(key), str(key)): count for key, count in histogram.items()}


def ratio(numerator: int, denominator: int) -> float:

if denominator == 0:
return 0.0

return numerator / denominator


def mib(byte_count: int) -> float:

return byte_count / (1024 * 1024)


def kib(byte_count: int) -> float:

return byte_count / 1024
25 changes: 0 additions & 25 deletions benchmark/suite/report.py
Original file line number Diff line number Diff line change
Expand Up @@ -213,34 +213,9 @@ def build_compact_report(

def memory_layout_rows(memory_layout: dict[str, Any]) -> list[dict[str, str]]:

instrument = memory_layout["instrument_sampling"]
optical = memory_layout["optical_accumulation"]
support = instrument["support_count_stats"]

return [
{
"area": "Instrument sampling plan",
"shape": (
f"{instrument['nominal_sample_count']} wavelengths, "
f"{instrument['kernel_count']} channel kernels, "
f"support count median {support['median']:.1f}, max {support['max']}"
),
"memory": (
f"current owned plan {instrument['current_owned_wavelength_sampling_bytes']} B; "
"estimated row payload saved "
f"{instrument['estimated_row_payload_saved_mib']:.2f} MiB"
),
},
{
"area": "Integration kernel scratch",
"shape": (
f"{instrument['estimated_sampling_worker_count']} sampling worker(s), "
f"{instrument['integration_kernel_scratch_bytes_per_worker']} B each"
),
"memory": (
f"{instrument['integration_kernel_scratch_bytes_at_worker_count']} B transient"
),
},
{
"area": "Collision-complex profile cache",
"shape": (
Expand Down
96 changes: 0 additions & 96 deletions python/zdisamar/bindings/handles.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,16 +18,12 @@
)
from ..output.tables import (
AtmosphericBudget,
InstrumentResponseTable,
O2LineContributions,
OxygenCollisionInducedAbsorptionDiagnosticTable,
)
from .loader import load_library
from .signatures import configure
from .structures import (
CAtmosphericBudget,
CDiagnosticReport,
CInstrumentResponse,
COptimalEstimationBatchRequest,
COptimalEstimationBatchResult,
COptimalEstimationControls,
Expand All @@ -37,8 +33,6 @@
COptimalEstimationResult,
COptimalEstimationScalarSpec,
CSpectrum,
O2LineContributionsRaw,
OxygenCollisionInducedAbsorptionDiagnosticsRaw,
)

_MAX_OPTIMAL_ESTIMATION_ITERATIONS = 1000
Expand All @@ -54,24 +48,6 @@ def contiguous_wavelengths(wavelengths_nm):
return double_array(wavelengths_nm, "wavelengths_nm")


def channel_mask(channels: tuple[str, ...]) -> int:
"""Translate requested spectral channels into model channel selection."""

masks = {"radiance": 1, "irradiance": 2}
mask = 0

for channel in channels:
try:
mask |= masks[channel]
except KeyError as exc:
raise ValueError(f"unsupported spectral channel: {channel}") from exc

if mask == 0:
raise ValueError("channels must not be empty")

return mask


def batch_run_status(value: int) -> str:
"""Translate native per-start batch status into a stable Python label."""

Expand Down Expand Up @@ -229,78 +205,6 @@ def atmospheric_budget(self, wavelengths_nm) -> AtmosphericBudget:

return AtmosphericBudget(self._copied_rows(raw, self._lib.zds_atmospheric_budget_free))

def o2_line_contributions(self, wavelengths_nm, max_rows: int = 50_000) -> O2LineContributions:
"""Return copied line-by-line evidence rows."""

wavelengths = contiguous_wavelengths(wavelengths_nm)

if max_rows <= 0:
raise ValueError("max_rows must be positive")

raw = O2LineContributionsRaw()
self._check(
self._lib.zds_o2_line_contributions(
self._ctx,
wavelengths,
len(wavelengths),
max_rows,
ctypes.byref(raw),
)
)
total_row_count = int(raw.total_row_count)
truncated = bool(raw.truncated)
rows = self._copied_rows(raw, self._lib.zds_o2_line_contributions_free)

return O2LineContributions(
rows,
total_row_count=total_row_count,
truncated=truncated,
)

def instrument_response_sampling(
self,
wavelengths_nm,
channels: tuple[str, ...] = ("radiance", "irradiance"),
) -> InstrumentResponseTable:
"""Return copied instrument response support rows."""

wavelengths = contiguous_wavelengths(wavelengths_nm)
raw = CInstrumentResponse()
self._check(
self._lib.zds_instrument_response_sampling(
self._ctx,
wavelengths,
len(wavelengths),
channel_mask(channels),
ctypes.byref(raw),
)
)

return InstrumentResponseTable(
self._copied_rows(raw, self._lib.zds_instrument_response_free)
)

def collision_induced_absorption(
self,
wavelengths_nm,
) -> OxygenCollisionInducedAbsorptionDiagnosticTable:
"""Return copied O2-O2 CIA rows on the atmospheric layer grid."""

wavelengths = contiguous_wavelengths(wavelengths_nm)
raw = OxygenCollisionInducedAbsorptionDiagnosticsRaw()
self._check(
self._lib.zds_o2_o2_cia_diagnostics(
self._ctx,
wavelengths,
len(wavelengths),
ctypes.byref(raw),
)
)

return OxygenCollisionInducedAbsorptionDiagnosticTable(
self._copied_rows(raw, self._lib.zds_o2_o2_cia_diagnostics_free)
)

def optimal_estimation(self, *, measurement, state_vector, controls):
"""Run native optimal estimation for the loaded scene."""

Expand Down
Loading