xlog_prob/mc/
mod.rs

1//! Approximate probabilistic inference via Monte Carlo sampling.
2//!
3//! This engine samples probabilistic facts / annotated disjunction decisions on the GPU and
4//! evaluates the deterministic core in each sampled world.
5//!
6//! For programs with non-monotone recursion (cycles through `not` and/or aggregates),
7//! Monte Carlo evaluation requires the user to opt into the approximate probabilistic
8//! engine. The deterministic evaluation uses a bounded, cycle-aware semantics:
9//!
10//! - If an SCC reaches a fixpoint under synchronous iteration, that fixpoint is used.
11//! - If the SCC enters a cycle, the interpretation is the intersection of all states in the cycle
12//!   (skeptical, invariant tuples only). This avoids parity/oscillation dependence on iteration
13//!   count while remaining fully deterministic and explicit.
14
15mod buffers;
16mod evidence;
17mod resident;
18mod results;
19
20pub use evidence::{EvidenceForcing, ForceabilityReason};
21pub use resident::{
22    compile_resident_plan, McNoHostStats, McResidentResult, ResidentPlan, ResidentRejectKind,
23    ResidentRejection,
24};
25
26use std::collections::BTreeMap;
27#[cfg(feature = "host-io")]
28use std::collections::{HashMap, HashSet};
29use std::sync::Arc;
30
31#[cfg(feature = "host-io")]
32use xlog_core::Schema;
33use xlog_core::{MemoryBudget, Result, XlogError};
34use xlog_cuda::memory::TrackedCudaSlice;
35use xlog_cuda::CudaKernelProvider;
36#[cfg(feature = "host-io")]
37use xlog_logic::ast::{BodyLiteral, Rule};
38use xlog_logic::ast::{Directives, Evidence, ProbMethod, ProbQuery, Program};
39
40use crate::exact::GpuConfig;
41#[cfg(feature = "host-io")]
42use crate::provenance::presentation_atom_from_canonical;
43#[cfg(feature = "host-io")]
44use crate::provenance::Value;
45use crate::provenance::{
46    atom_key_from_ground_atom, canonicalize_probabilistic_program, GroundAtom,
47};
48
49/// Sampling method for Monte Carlo inference.
50#[derive(Debug, Clone, Copy, PartialEq, Eq)]
51pub enum McSamplingMethod {
52    /// Sample from prior, discard worlds where evidence is not satisfied.
53    Rejection,
54    /// Force evidence variables in the sampler; every sample counts.
55    EvidenceClamping,
56}
57
58impl McSamplingMethod {
59    /// Returns the stable machine-readable sampling-method name.
60    pub fn as_str(self) -> &'static str {
61        match self {
62            McSamplingMethod::Rejection => "rejection",
63            McSamplingMethod::EvidenceClamping => "evidence_clamping",
64        }
65    }
66}
67
68impl From<ProbMethod> for McSamplingMethod {
69    fn from(value: ProbMethod) -> Self {
70        match value {
71            ProbMethod::Rejection => Self::Rejection,
72            ProbMethod::EvidenceClamping => Self::EvidenceClamping,
73        }
74    }
75}
76
77/// Strategy for counting evidence-satisfied samples in the MC loop.
78///
79/// In `QueriesOnly` mode (used with evidence clamping), evidence is
80/// guaranteed to hold in every sample, so we skip the truth-kernel's
81/// evidence check and evidence-side buffer allocations.
82#[derive(Debug, Clone, Copy, PartialEq, Eq)]
83pub enum McCountStrategy {
84    /// Full path: check both queries and evidence each sample.
85    QueriesAndEvidence,
86    /// Clamped path: evidence is always satisfied; only accumulate query flags.
87    QueriesOnly,
88}
89
90impl McCountStrategy {
91    /// Derive the count strategy from the chosen sampling method.
92    pub fn from_method(method: McSamplingMethod) -> Self {
93        match method {
94            McSamplingMethod::Rejection => Self::QueriesAndEvidence,
95            McSamplingMethod::EvidenceClamping => Self::QueriesOnly,
96        }
97    }
98}
99
100/// Breakdown of time spent in each phase of MC evaluation.
101/// Gate with `XLOG_MC_PROFILE=1` to print at the end of evaluation.
102#[derive(Debug, Clone, Copy, Default, PartialEq, Eq)]
103pub struct McTimingBreakdown {
104    /// Microseconds spent drawing Bernoulli samples.
105    pub sampler_us: u64,
106    /// Microseconds spent resetting per-sample state.
107    pub sample_reset_us: u64,
108    /// Microseconds spent materializing sampled extensional facts.
109    pub sample_build_us: u64,
110    /// Microseconds spent evaluating deterministic rules.
111    pub eval_us: u64,
112    /// Microseconds spent accumulating evidence and query counts.
113    pub count_us: u64,
114}
115
116impl McTimingBreakdown {
117    /// Returns the saturating sum of all measured phases.
118    pub fn total_us(&self) -> u64 {
119        self.sampler_us
120            .saturating_add(self.sample_reset_us)
121            .saturating_add(self.sample_build_us)
122            .saturating_add(self.eval_us)
123            .saturating_add(self.count_us)
124    }
125}
126
127/// Bounded semantics for non-monotone SCC evaluation inside MC sampling.
128pub const NONMONOTONE_SEMANTICS: &str = "Synchronous iteration per SCC; if a fixpoint is reached, use it; if a cycle is detected, use the intersection of all states in the cycle (skeptical tuples only); if the iteration budget is exceeded, use the intersection across all visited states (conservative).";
129
130#[derive(Debug, Clone)]
131#[non_exhaustive]
132/// Configuration for Monte Carlo probabilistic inference.
133///
134/// Use [`McEvalConfig::default()`] as a starting point and then update the
135/// individual fields you need.
136pub struct McEvalConfig {
137    /// Number of Monte Carlo samples.
138    pub samples: usize,
139    /// RNG seed (deterministic).
140    pub seed: u64,
141    /// Two-sided confidence level in (0,1) (e.g., 0.95).
142    pub confidence: f64,
143    /// Maximum SCC iteration steps for non-monotone cycle detection.
144    pub max_nonmonotone_iterations: usize,
145    /// Sampling method override. `None` = auto-select (EvidenceClamping when forceable, Rejection otherwise).
146    pub sampling_method: Option<McSamplingMethod>,
147    /// Allow the host CPU oracle ([`McProgram::evaluate_cpu`]) when the
148    /// resident GPU engine rejects the program. Default `false`: a rejected
149    /// program fails closed with the typed rejection instead of silently
150    /// running on the CPU. When set, the oracle result is labeled
151    /// [`McEngine::CpuOracle`] and must never serve as GPU-native evidence.
152    pub allow_cpu_oracle_fallback: bool,
153}
154
155/// Explicit caller overrides layered on top of source directives.
156///
157/// `None` preserves the corresponding source directive or language default.
158#[derive(Debug, Clone, Default)]
159pub struct McEvalOverrides {
160    /// Replacement sample count.
161    pub samples: Option<usize>,
162    /// Replacement deterministic seed.
163    pub seed: Option<u64>,
164    /// Replacement two-sided confidence level.
165    pub confidence: Option<f64>,
166    /// Replacement non-monotone iteration limit.
167    pub max_nonmonotone_iterations: Option<usize>,
168    /// Replacement sampling method.
169    pub sampling_method: Option<McSamplingMethod>,
170    /// Replacement explicit CPU-oracle opt-in.
171    pub allow_cpu_oracle_fallback: Option<bool>,
172}
173
174impl Default for McEvalConfig {
175    fn default() -> Self {
176        Self {
177            samples: 10000,
178            seed: 0,
179            confidence: 0.95,
180            max_nonmonotone_iterations: 1024,
181            sampling_method: None,
182            allow_cpu_oracle_fallback: false,
183        }
184    }
185}
186
187impl McEvalConfig {
188    /// Resolves Monte Carlo settings from program directives and validates them.
189    pub fn from_directives(directives: &Directives) -> Result<Self> {
190        let mut cfg = Self::default();
191        if let Some(samples) = directives.prob_samples {
192            cfg.samples = samples;
193        }
194        if let Some(seed) = directives.prob_seed {
195            cfg.seed = seed;
196        }
197        if let Some(confidence) = directives.prob_confidence {
198            cfg.confidence = confidence;
199        }
200        if let Some(iterations) = directives.prob_max_nonmonotone_iterations {
201            cfg.max_nonmonotone_iterations = iterations;
202        }
203        cfg.sampling_method = directives.prob_method.map(McSamplingMethod::from);
204        cfg.validate()?;
205        Ok(cfg)
206    }
207
208    /// Apply only explicitly supplied caller values, then validate the resolved config.
209    pub fn apply_overrides(&mut self, overrides: McEvalOverrides) -> Result<()> {
210        if let Some(samples) = overrides.samples {
211            self.samples = samples;
212        }
213        if let Some(seed) = overrides.seed {
214            self.seed = seed;
215        }
216        if let Some(confidence) = overrides.confidence {
217            self.confidence = confidence;
218        }
219        if let Some(iterations) = overrides.max_nonmonotone_iterations {
220            self.max_nonmonotone_iterations = iterations;
221        }
222        if let Some(method) = overrides.sampling_method {
223            self.sampling_method = Some(method);
224        }
225        if let Some(allow) = overrides.allow_cpu_oracle_fallback {
226            self.allow_cpu_oracle_fallback = allow;
227        }
228        self.validate()
229    }
230
231    /// Validates numeric ranges and required positive limits.
232    pub fn validate(&self) -> Result<()> {
233        if self.samples == 0 {
234            return Err(XlogError::Compilation(
235                "MC inference requires samples > 0".to_string(),
236            ));
237        }
238        if !(0.0 < self.confidence && self.confidence < 1.0) || self.confidence.is_nan() {
239            return Err(XlogError::Compilation(format!(
240                "MC inference requires 0 < confidence < 1, got {}",
241                self.confidence
242            )));
243        }
244        if self.max_nonmonotone_iterations == 0 {
245            return Err(XlogError::Compilation(
246                "MC inference requires max_nonmonotone_iterations > 0".to_string(),
247            ));
248        }
249        Ok(())
250    }
251}
252
253#[derive(Debug, Clone)]
254/// Monte Carlo probability estimate and confidence interval for one query.
255pub struct McQueryEstimate {
256    /// Ground query atom.
257    pub atom: GroundAtom,
258    /// Estimated conditional probability.
259    pub prob: f64,
260    /// Natural logarithm of the estimated conditional probability.
261    pub log_prob: f64,
262    /// Estimated standard error.
263    pub stderr: f64,
264    /// Lower endpoint of the configured confidence interval.
265    pub ci_low: f64,
266    /// Upper endpoint of the configured confidence interval.
267    pub ci_high: f64,
268}
269
270/// Which engine produced an [`McResult`].
271#[derive(Debug, Clone, Copy, PartialEq, Eq)]
272pub enum McEngine {
273    /// GPU-resident megakernel engine — the production MC path.
274    GpuResident,
275    /// Host CPU oracle ([`McProgram::evaluate_cpu`]) — explicit opt-in only
276    /// (see [`McEvalConfig::allow_cpu_oracle_fallback`]); never valid as
277    /// GPU-native or zero-host evidence.
278    CpuOracle,
279}
280
281impl McEngine {
282    /// Stable string form for result metadata surfaces (CLI JSON, pyxlog).
283    pub fn as_str(&self) -> &'static str {
284        match self {
285            McEngine::GpuResident => "gpu-resident",
286            McEngine::CpuOracle => "cpu-oracle",
287        }
288    }
289}
290
291#[derive(Debug, Clone)]
292/// Materialized Monte Carlo evaluation result.
293pub struct McResult {
294    /// Number of sampled worlds.
295    pub total_samples: usize,
296    /// Number of sampled worlds satisfying all evidence.
297    pub evidence_samples: usize,
298    /// Random seed used for sampling.
299    pub seed: u64,
300    /// Confidence level used for intervals.
301    pub confidence: f64,
302    /// Per-query conditional estimates.
303    pub query_estimates: Vec<McQueryEstimate>,
304    /// Number of non-monotone strongly connected components evaluated.
305    pub nonmonotone_sccs: usize,
306    /// Number of detected non-monotone evaluation cycles.
307    pub nonmonotone_cycles: usize,
308    /// Number of non-monotone components that reached the iteration limit.
309    pub nonmonotone_iteration_limit_hits: usize,
310    /// Evidence-sampling method used for the run.
311    pub sampling_method: McSamplingMethod,
312    /// Engine that produced this result; CPU-oracle results are reachable
313    /// only through explicit opt-in and are labeled so downstream consumers
314    /// can never mistake them for GPU-resident evidence.
315    pub engine: McEngine,
316}
317
318/// **Legacy back-compat surface** — tracked (data-plane) host<->device transfer
319/// deltas measured around the MC measured region.
320///
321/// This struct dates from the predecessor `a894aab4` engine, which removed only
322/// *tracked* data-plane transfers from a still-host-orchestrated loop. The
323/// current resident megakernel engine has a *stronger* property — **no host
324/// interaction at all** in the measured region — so the authoritative contract
325/// now lives in [`McNoHostStats`] (`McResidentResult::no_host`), which also
326/// counts untracked metadata reads, host fixpoint iterations, and in-region
327/// device allocations. `McDeviceResult` retains this field for API
328/// back-compatibility; for the resident engine its tracked-call fields are zero.
329#[derive(Debug, Clone, Copy, Default, PartialEq, Eq)]
330pub struct McHotLoopTransfers {
331    /// Tracked host-to-device calls in the measured region.
332    pub htod_calls: u64,
333    /// Tracked device-to-host calls in the measured region.
334    pub dtoh_calls: u64,
335    /// Tracked host-to-device bytes in the measured region.
336    pub htod_bytes: u64,
337    /// Tracked device-to-host bytes in the measured region.
338    pub dtoh_bytes: u64,
339}
340
341impl McHotLoopTransfers {
342    /// True when no tracked host/device transfer occurred in the measured region.
343    pub fn is_zero(&self) -> bool {
344        self.htod_calls == 0 && self.dtoh_calls == 0
345    }
346}
347
348/// Device-resident Monte Carlo result counts.
349pub struct McDeviceResult {
350    /// Device count of evidence-satisfying worlds that also satisfy each query.
351    pub query_counts: TrackedCudaSlice<u32>,
352    /// Device count of worlds satisfying all evidence.
353    pub evidence_count: TrackedCudaSlice<u32>,
354    /// Number of sampled worlds.
355    pub total_samples: usize,
356    /// Random seed used for sampling.
357    pub seed: u64,
358    /// Confidence level requested for result materialization.
359    pub confidence: f64,
360    /// Number of non-monotone strongly connected components evaluated.
361    pub nonmonotone_sccs: usize,
362    /// Number of detected non-monotone evaluation cycles.
363    pub nonmonotone_cycles: usize,
364    /// Number of non-monotone components that reached the iteration limit.
365    pub nonmonotone_iteration_limit_hits: usize,
366    /// Evidence-sampling method used for the run.
367    pub sampling_method: McSamplingMethod,
368    /// Legacy back-compat field: tracked transfers measured around the resident
369    /// engine's measured region (zero). The authoritative no-host contract is
370    /// [`McNoHostStats`] on [`McResidentResult`]; see [`McHotLoopTransfers`].
371    pub hot_loop_transfers: McHotLoopTransfers,
372    /// Authoritative no-host counters for the resident engine measured region.
373    pub no_host: McNoHostStats,
374}
375
376#[derive(Debug, Clone)]
377pub(super) struct ProbFactSpec {
378    pub(super) var_idx: usize,
379    pub(super) atom: GroundAtom,
380}
381
382#[derive(Debug, Clone)]
383pub(super) struct AdSpec {
384    pub(super) decision_vars: Vec<usize>,
385    pub(super) choices: Vec<GroundAtom>,
386    #[cfg(feature = "host-io")]
387    pub(super) has_none: bool,
388}
389
390#[cfg(feature = "host-io")]
391#[derive(Debug, Clone, Default, PartialEq, Eq)]
392pub(super) struct Relation {
393    pub(super) tuples: HashSet<Vec<Value>>,
394}
395
396#[cfg(feature = "host-io")]
397impl Relation {
398    pub(super) fn insert_tuple(&mut self, tuple: Vec<Value>) {
399        self.tuples.insert(tuple);
400    }
401
402    pub(super) fn contains(&self, tuple: &[Value]) -> bool {
403        self.tuples.contains(tuple)
404    }
405
406    pub(super) fn is_empty(&self) -> bool {
407        self.tuples.is_empty()
408    }
409}
410
411#[cfg(feature = "host-io")]
412#[derive(Debug, Clone)]
413pub(super) enum SccKind {
414    MonotoneNonRecursive,
415    MonotoneRecursive,
416    NonMonotone,
417}
418
419#[cfg(feature = "host-io")]
420#[derive(Debug, Clone)]
421pub(super) struct SccPlan {
422    pub(super) predicates: Vec<String>,
423    pub(super) rules: Vec<Rule>,
424    pub(super) kind: SccKind,
425}
426
427#[cfg(feature = "host-io")]
428#[derive(Debug, Clone, Default)]
429pub(super) struct EvalStats {
430    pub(super) nonmonotone_sccs: usize,
431    pub(super) nonmonotone_cycles: usize,
432    pub(super) nonmonotone_iteration_limit_hits: usize,
433}
434
435#[derive(Clone)]
436/// Compiled Monte Carlo program with canonical probabilistic variables and GPU settings.
437pub struct McProgram {
438    pub(super) gpu_config: GpuConfig,
439    pub(super) program: Program,
440    #[cfg(feature = "host-io")]
441    pub(super) base_store: HashMap<String, Relation>,
442    #[cfg(feature = "host-io")]
443    pub(super) scc_plans: Vec<SccPlan>,
444    #[cfg(feature = "host-io")]
445    pub(super) arithmetic_schemas: HashMap<String, Schema>,
446    pub(super) queries: Vec<GroundAtom>,
447    #[cfg(feature = "host-io")]
448    query_presentations: Vec<GroundAtom>,
449    pub(super) evidence: Vec<(GroundAtom, bool)>,
450    pub(super) bernoulli_probs: Vec<f32>,
451    pub(super) prob_facts: Vec<ProbFactSpec>,
452    pub(super) annotated_disjunctions: Vec<AdSpec>,
453}
454
455impl McProgram {
456    /// Parses source and compiles a Monte Carlo program with default GPU settings.
457    pub fn compile_source(source: &str) -> Result<Self> {
458        let program = xlog_logic::parse_program(source)?;
459        Self::compile_from_program(&program, GpuConfig::default())
460    }
461
462    /// Parses source and compiles a Monte Carlo program with explicit GPU settings.
463    pub fn compile_source_with_gpu(source: &str, config: GpuConfig) -> Result<Self> {
464        let program = xlog_logic::parse_program(source)?;
465        Self::compile_from_program(&program, config)
466    }
467
468    /// Compile an already parsed program with the requested GPU configuration.
469    ///
470    /// Imports must already be resolved and merged. This method does not load
471    /// unresolved `use` declarations from the filesystem.
472    pub fn compile_from_program(program: &Program, config: GpuConfig) -> Result<Self> {
473        let mut compiled = Self::compile_program(program)?;
474        compiled.gpu_config = config;
475        Ok(compiled)
476    }
477
478    /// Returns the number of independent Bernoulli variables sampled per world.
479    pub fn num_vars(&self) -> usize {
480        self.bernoulli_probs.len()
481    }
482
483    /// Host-facing MC evaluation: runs the resident megakernel engine
484    /// ([`Self::evaluate_gpu_device_with_provider`]) and then **materializes the
485    /// result on the host** by downloading the final query/evidence counts
486    /// *after* the measured region. The download is a host-result
487    /// materialization, not part of the measured region — the no-host property
488    /// belongs to the resident engine, not to this convenience wrapper. Use
489    /// [`Self::evaluate_gpu_device`] when you want device-resident counts with no
490    /// host download at all.
491    ///
492    /// **Fail-closed contract:** if the resident engine rejects the program
493    /// (negation, aggregates, unbounded terms, ...), this returns the typed
494    /// rejection error. The CPU oracle is reachable only via the explicit
495    /// [`McEvalConfig::allow_cpu_oracle_fallback`] opt-in and its result is
496    /// labeled [`McEngine::CpuOracle`].
497    #[cfg(feature = "host-io")]
498    pub fn evaluate(&self, cfg: McEvalConfig) -> Result<McResult> {
499        let provider = Arc::new(self.provider()?);
500        let cfg_clone = cfg.clone();
501        let device_result =
502            match self.evaluate_gpu_device_with_provider(cfg_clone, provider.clone()) {
503                Ok(result) => result,
504                Err(XlogError::Compilation(message))
505                    if message.starts_with("resident MC engine rejected program") =>
506                {
507                    // Fail closed by default: the CPU oracle is an explicit,
508                    // labeled opt-in (`allow_cpu_oracle_fallback`), never a
509                    // silent substitute for the resident GPU engine.
510                    if cfg.allow_cpu_oracle_fallback {
511                        return self.evaluate_cpu(cfg);
512                    }
513                    return Err(XlogError::Compilation(format!(
514                        "{message}; MC fail-closed: set \
515                         McEvalConfig::allow_cpu_oracle_fallback to explicitly \
516                         run the labeled CPU oracle instead"
517                    )));
518                }
519                Err(err) => return Err(err),
520            };
521
522        let mut host_counts = vec![0u32; device_result.query_counts.len()];
523        if !host_counts.is_empty() {
524            provider
525                .device()
526                .inner()
527                .dtoh_sync_copy_into(&device_result.query_counts, &mut host_counts)
528                .map_err(|e| {
529                    XlogError::Kernel(format!("Failed to download MC query counts: {}", e))
530                })?;
531        }
532
533        let mut host_evidence = [0u32];
534        provider
535            .device()
536            .inner()
537            .dtoh_sync_copy_into(&device_result.evidence_count, &mut host_evidence)
538            .map_err(|e| {
539                XlogError::Kernel(format!("Failed to download MC evidence count: {}", e))
540            })?;
541
542        let evidence_samples = if self.evidence.is_empty() {
543            cfg.samples
544        } else {
545            host_evidence[0] as usize
546        };
547
548        if device_result.sampling_method != McSamplingMethod::EvidenceClamping
549            && !self.evidence.is_empty()
550            && evidence_samples == 0
551        {
552            return Err(XlogError::Execution(format!(
553                "MC inference error: evidence was never satisfied across {} samples (seed={})",
554                cfg.samples, cfg.seed
555            )));
556        }
557
558        let z = results::normal_quantile(0.5 + cfg.confidence / 2.0);
559        let mut query_estimates: Vec<McQueryEstimate> = Vec::with_capacity(self.queries.len());
560        for (i, atom) in self.query_presentations.iter().enumerate() {
561            let k = host_counts.get(i).copied().unwrap_or(0) as usize;
562            let (p, stderr, ci_low, ci_high) = results::binomial_estimate(k, evidence_samples, z);
563            let log_prob = if p == 0.0 { f64::NEG_INFINITY } else { p.ln() };
564
565            query_estimates.push(McQueryEstimate {
566                atom: atom.clone(),
567                prob: p,
568                log_prob,
569                stderr,
570                ci_low,
571                ci_high,
572            });
573        }
574
575        Ok(McResult {
576            total_samples: cfg.samples,
577            evidence_samples,
578            seed: cfg.seed,
579            confidence: cfg.confidence,
580            query_estimates,
581            nonmonotone_sccs: device_result.nonmonotone_sccs,
582            nonmonotone_cycles: device_result.nonmonotone_cycles,
583            engine: McEngine::GpuResident,
584            nonmonotone_iteration_limit_hits: device_result.nonmonotone_iteration_limit_hits,
585            sampling_method: device_result.sampling_method,
586        })
587    }
588
589    /// CPU **oracle / debug** MC path. Downloads the full sampled-bit matrix to
590    /// the host and evaluates every sampled world on a host relation store.
591    ///
592    /// This is intentionally *not* GPU-native: it performs a large DtoH of the
593    /// sample matrix and runs the deterministic core on the CPU. It exists solely
594    /// as a deterministic, seed-matched oracle for validating the GPU-native
595    /// device counts (the GPU sampler is shared, so for the same program/seed the
596    /// two paths see identical samples). It must **never** be used as zero-host /
597    /// GPU-native release evidence, and the acceptance matrix excludes it and the
598    /// tests that call it (`tests/gpu_mc_vs_cpu.rs`, `tests/mc.rs`).
599    #[cfg(feature = "host-io")]
600    pub fn evaluate_cpu(&self, cfg: McEvalConfig) -> Result<McResult> {
601        if cfg.samples == 0 {
602            return Err(XlogError::Execution(
603                "MC inference requires samples > 0".to_string(),
604            ));
605        }
606        if !(0.0 < cfg.confidence && cfg.confidence < 1.0) || cfg.confidence.is_nan() {
607            return Err(XlogError::Execution(format!(
608                "MC inference requires 0 < confidence < 1, got {}",
609                cfg.confidence
610            )));
611        }
612        if cfg.max_nonmonotone_iterations == 0 {
613            return Err(XlogError::Execution(
614                "MC inference requires max_nonmonotone_iterations > 0".to_string(),
615            ));
616        }
617
618        let (method, forcing) = self.resolve_sampling_method(cfg.sampling_method)?;
619        let is_clamped = method == McSamplingMethod::EvidenceClamping;
620
621        let mut n_evidence: usize = 0;
622        let mut n_query_true: Vec<usize> = vec![0; self.queries.len()];
623        let mut stats = EvalStats::default();
624
625        let num_vars = self.bernoulli_probs.len();
626        let samples_matrix: Vec<u8> = if num_vars == 0 {
627            Vec::new()
628        } else {
629            let total = num_vars
630                .checked_mul(cfg.samples)
631                .ok_or_else(|| XlogError::Execution("MC sample matrix overflow".to_string()))?;
632            let provider = Arc::new(self.provider()?);
633
634            // Allocate force arrays: upload actual forcing data in clamped mode, zero-fill otherwise
635            let mut d_force_mask = provider.memory().alloc::<u8>(num_vars.max(1))?;
636            let mut d_forced_value = provider.memory().alloc::<u8>(num_vars.max(1))?;
637            if is_clamped {
638                provider
639                    .htod_sync_copy_into_tracked(&forcing.force_mask, &mut d_force_mask)
640                    .map_err(|e| {
641                        XlogError::Kernel(format!("Failed to upload force_mask: {}", e))
642                    })?;
643                provider
644                    .htod_sync_copy_into_tracked(&forcing.forced_value, &mut d_forced_value)
645                    .map_err(|e| {
646                        XlogError::Kernel(format!("Failed to upload forced_value: {}", e))
647                    })?;
648            } else {
649                provider
650                    .device()
651                    .inner()
652                    .memset_zeros(&mut d_force_mask)
653                    .map_err(|e| XlogError::Kernel(format!("Failed to zero force_mask: {}", e)))?;
654                provider
655                    .device()
656                    .inner()
657                    .memset_zeros(&mut d_forced_value)
658                    .map_err(|e| {
659                        XlogError::Kernel(format!("Failed to zero forced_value: {}", e))
660                    })?;
661            }
662
663            let samples_device = provider.sample_bernoulli_matrix_device(
664                &self.bernoulli_probs,
665                cfg.samples,
666                cfg.seed,
667                &d_force_mask.slice(..),
668                &d_forced_value.slice(..),
669            )?;
670            let mut host = vec![0u8; total];
671            if !host.is_empty() {
672                provider
673                    .device()
674                    .inner()
675                    .dtoh_sync_copy_into(&samples_device, &mut host)
676                    .map_err(|e| {
677                        XlogError::Kernel(format!("Failed to download MC samples: {}", e))
678                    })?;
679            }
680            host
681        };
682
683        for sample_idx in 0..cfg.samples {
684            let sample_bits = if num_vars == 0 {
685                &[][..]
686            } else {
687                let start = sample_idx * num_vars;
688                let end = start + num_vars;
689                &samples_matrix[start..end]
690            };
691
692            let mut store = self.base_store.clone();
693            self.apply_sample_facts(&mut store, sample_bits)?;
694
695            let sample_stats = results::evaluate_program_inplace(
696                &self.scc_plans,
697                &mut store,
698                cfg.max_nonmonotone_iterations,
699                &self.arithmetic_schemas,
700            )?;
701            stats.nonmonotone_sccs += sample_stats.nonmonotone_sccs;
702            stats.nonmonotone_cycles += sample_stats.nonmonotone_cycles;
703            stats.nonmonotone_iteration_limit_hits += sample_stats.nonmonotone_iteration_limit_hits;
704
705            // In clamped mode, skip evidence check — all samples count
706            if !is_clamped && !results::evidence_satisfied(&store, &self.evidence) {
707                continue;
708            }
709
710            n_evidence += 1;
711            for (i, q) in self.queries.iter().enumerate() {
712                if results::atom_holds(&store, q) {
713                    n_query_true[i] += 1;
714                }
715            }
716        }
717
718        if !is_clamped && !self.evidence.is_empty() && n_evidence == 0 {
719            return Err(XlogError::Execution(format!(
720                "MC inference error: evidence was never satisfied across {} samples (seed={})",
721                cfg.samples, cfg.seed
722            )));
723        }
724
725        // If there is no evidence (or clamped mode), treat all samples as evidence-satisfying.
726        let denom = if self.evidence.is_empty() || is_clamped {
727            cfg.samples
728        } else {
729            n_evidence
730        };
731
732        let z = results::normal_quantile(0.5 + cfg.confidence / 2.0);
733
734        let mut query_estimates: Vec<McQueryEstimate> = Vec::with_capacity(self.queries.len());
735        for (i, atom) in self.query_presentations.iter().enumerate() {
736            let k = n_query_true[i];
737            let (p, stderr, ci_low, ci_high) = results::binomial_estimate(k, denom, z);
738            let log_prob = if p == 0.0 { f64::NEG_INFINITY } else { p.ln() };
739
740            query_estimates.push(McQueryEstimate {
741                atom: atom.clone(),
742                prob: p,
743                log_prob,
744                stderr,
745                ci_low,
746                ci_high,
747            });
748        }
749
750        Ok(McResult {
751            total_samples: cfg.samples,
752            evidence_samples: denom,
753            seed: cfg.seed,
754            confidence: cfg.confidence,
755            query_estimates,
756            nonmonotone_sccs: stats.nonmonotone_sccs,
757            nonmonotone_cycles: stats.nonmonotone_cycles,
758            nonmonotone_iteration_limit_hits: stats.nonmonotone_iteration_limit_hits,
759            sampling_method: method,
760            engine: McEngine::CpuOracle,
761        })
762    }
763
764    /// Alias for [`Self::evaluate`]: GPU device evaluation followed by host-result
765    /// materialization (final-count download after the measured region). The
766    /// `_gpu` suffix denotes that the *compute* runs on the GPU — it does **not**
767    /// imply a zero-host result, since it returns a host [`McResult`]. For the
768    /// device-resident, no-host-download API use [`Self::evaluate_gpu_device`].
769    #[cfg(feature = "host-io")]
770    pub fn evaluate_gpu(&self, cfg: McEvalConfig) -> Result<McResult> {
771        self.evaluate(cfg)
772    }
773
774    /// GPU-native device-resident MC evaluation via the resident megakernel
775    /// engine ([`resident`]). Returns [`McDeviceResult`] with counts left on the
776    /// device (no host download). The engine evaluates all worlds in a single
777    /// launch with **no host interaction in the measured region** (no host
778    /// sample loop, no per-sample/per-operator host launches or allocations, no
779    /// tracked transfers, and no untracked metadata reads); see
780    /// [`McResidentResult::no_host`] / [`McNoHostStats::is_no_host`] for the full
781    /// measured contract.
782    pub fn evaluate_gpu_device(&self, cfg: McEvalConfig) -> Result<McDeviceResult> {
783        let provider = Arc::new(self.provider()?);
784        self.evaluate_gpu_device_with_provider(cfg, provider)
785    }
786
787    /// GPU-native device-resident MC evaluation using a caller-supplied provider
788    /// (enables provider/buffer reuse across calls). Static setup (arena
789    /// allocation, plan upload, sampling) happens before the measured region;
790    /// the measured region is a single resident-engine launch with **zero host
791    /// interaction** — no host loop over samples or fixpoint iterations, no
792    /// per-sample/per-operator host launches or device allocations, no tracked
793    /// HtoD/DtoH transfers, and no untracked metadata reads. The full contract is
794    /// measured by [`McNoHostStats`] (`McResidentResult::no_host`) and gated by
795    /// `tests/mc_resident.rs`. Counts remain device-resident; the caller decides
796    /// whether/when to download them.
797    pub fn evaluate_gpu_device_with_provider(
798        &self,
799        cfg: McEvalConfig,
800        provider: Arc<CudaKernelProvider>,
801    ) -> Result<McDeviceResult> {
802        // The GPU-resident megakernel engine is the sole MC execution path: there
803        // is no host-orchestrated per-sample fallback. It evaluates ALL worlds in
804        // a single launch with zero host interaction in the measured region, then
805        // returns device-resident counts. Programs outside the supported bounded
806        // fragment fail closed (typed `ResidentRejection`).
807        let r = self.evaluate_resident_with_provider(cfg, provider)?;
808        let no_host = r.no_host;
809        Ok(McDeviceResult {
810            query_counts: r.query_counts,
811            evidence_count: r.evidence_count,
812            total_samples: r.total_samples,
813            seed: r.seed,
814            confidence: r.confidence,
815            // The resident engine accepts only bounded positive Datalog and
816            // device-evaluates its fixpoint; legacy non-monotone SCC bookkeeping
817            // is not part of this engine's reported state.
818            nonmonotone_sccs: 0,
819            nonmonotone_cycles: 0,
820            nonmonotone_iteration_limit_hits: 0,
821            sampling_method: r.sampling_method,
822            // Back-compat surface: `hot_loop_transfers` reports the resident
823            // engine's tracked transfers (zero). Richer no-host evidence (untracked
824            // reads, fixpoint/alloc counts) lives in `McResidentResult::no_host`.
825            hot_loop_transfers: McHotLoopTransfers {
826                htod_calls: no_host.tracked_htod_calls,
827                dtoh_calls: no_host.tracked_dtoh_calls,
828                htod_bytes: 0,
829                dtoh_bytes: 0,
830            },
831            no_host,
832        })
833    }
834
835    fn compile_program(program: &Program) -> Result<Self> {
836        #[cfg(feature = "host-io")]
837        let source_program = program;
838        let arithmetic_schemas = crate::provenance::arithmetic_schemas(program)?;
839        let canonical_program = canonicalize_probabilistic_program(program, &arithmetic_schemas)?;
840        let program = &canonical_program;
841
842        let mut queries: Vec<GroundAtom> = Vec::new();
843        for ProbQuery { atom } in &program.prob_queries {
844            queries.push(atom_key_from_ground_atom(atom)?);
845        }
846        #[cfg(feature = "host-io")]
847        let mut query_presentations: Vec<GroundAtom> = Vec::new();
848        #[cfg(feature = "host-io")]
849        {
850            for (ProbQuery { atom }, ProbQuery { atom: source_atom }) in program
851                .prob_queries
852                .iter()
853                .zip(&source_program.prob_queries)
854            {
855                query_presentations.push(presentation_atom_from_canonical(
856                    source_atom,
857                    atom,
858                    &arithmetic_schemas,
859                )?);
860            }
861        }
862
863        let mut evidence: Vec<(GroundAtom, bool)> = Vec::new();
864        let mut evidence_values: BTreeMap<GroundAtom, bool> = BTreeMap::new();
865        for Evidence { atom, value } in &program.evidence {
866            let atom = atom_key_from_ground_atom(atom)?;
867            match evidence_values.entry(atom.clone()) {
868                std::collections::btree_map::Entry::Occupied(previous)
869                    if previous.get() != value =>
870                {
871                    return Err(XlogError::Execution(format!(
872                        "Monte Carlo inference error: conflicting evidence for {}",
873                        atom.predicate
874                    )));
875                }
876                std::collections::btree_map::Entry::Occupied(_) => continue,
877                std::collections::btree_map::Entry::Vacant(entry) => {
878                    entry.insert(*value);
879                    evidence.push((atom, *value));
880                }
881            }
882        }
883
884        let mut prob_facts = program.prob_facts.clone();
885        buffers::extend_prob_facts_with_coin(program, &mut prob_facts)?;
886        let (bernoulli_probs, prob_facts, annotated_disjunctions) =
887            buffers::compile_sampling_plan(&prob_facts, &program.annotated_disjunctions)?;
888
889        #[cfg(feature = "host-io")]
890        let mut base_store: HashMap<String, Relation> = HashMap::new();
891
892        // Deterministic facts.
893        #[cfg(feature = "host-io")]
894        {
895            for fact in program.facts() {
896                let atom = atom_key_from_ground_atom(&fact.head)?;
897                base_store
898                    .entry(atom.predicate.clone())
899                    .or_default()
900                    .insert_tuple(atom.args);
901            }
902        }
903
904        // Ensure relations exist for all referenced predicates so evaluation treats missing as empty,
905        // but never errors due to an unknown predicate (CPU eval path only).
906        #[cfg(feature = "host-io")]
907        {
908            let mut referenced: HashSet<String> = HashSet::new();
909            for rule in &program.rules {
910                referenced.insert(rule.head.predicate.clone());
911                for lit in &rule.body {
912                    match lit {
913                        BodyLiteral::Positive(a) | BodyLiteral::Negated(a) => {
914                            referenced.insert(a.predicate.clone());
915                        }
916                        BodyLiteral::Epistemic(lit) => {
917                            referenced.insert(lit.atom.predicate.clone());
918                        }
919                        BodyLiteral::Comparison(_)
920                        | BodyLiteral::IsExpr(_)
921                        | BodyLiteral::Univ(_) => {}
922                    }
923                }
924            }
925            for pf in &program.prob_facts {
926                referenced.insert(pf.atom.predicate.clone());
927            }
928            for ad in &program.annotated_disjunctions {
929                for pf in &ad.choices {
930                    referenced.insert(pf.atom.predicate.clone());
931                }
932            }
933            for q in &queries {
934                referenced.insert(q.predicate.clone());
935            }
936            for (e, _) in &evidence {
937                referenced.insert(e.predicate.clone());
938            }
939            for pred in referenced {
940                base_store.entry(pred).or_default();
941            }
942        }
943
944        #[cfg(feature = "host-io")]
945        let scc_plans = results::build_scc_plans(program)?;
946        Ok(Self {
947            gpu_config: GpuConfig::default(),
948            program: program.clone(),
949            #[cfg(feature = "host-io")]
950            base_store,
951            #[cfg(feature = "host-io")]
952            scc_plans,
953            #[cfg(feature = "host-io")]
954            arithmetic_schemas,
955            queries,
956            #[cfg(feature = "host-io")]
957            query_presentations,
958            evidence,
959            bernoulli_probs,
960            prob_facts,
961            annotated_disjunctions,
962        })
963    }
964
965    fn provider(&self) -> Result<CudaKernelProvider> {
966        if self.gpu_config.memory_bytes == 0 {
967            return Err(XlogError::Kernel(
968                "GPU memory budget must be non-zero".to_string(),
969            ));
970        }
971
972        xlog_cuda::CudaProviderBuilder::new(
973            self.gpu_config.device_ordinal,
974            MemoryBudget::with_limit(self.gpu_config.memory_bytes),
975        )
976        .build()
977    }
978}
xlog_prob/mc/mod.rs

xlog_prob/mc/
mod.rs
