GPU-beschleunigte Unsicherheit Quantifizierung: Skalierung von Monte-Carlo- und PCE-Methoden auf der GPU

Key Takeaways

• Die GPU-Beschleunigung bietet eine 10–1000-fache Beschleunigung für die Batch-Monte-Carlo-Sampling im Vergleich zu nur CPU-Ansätzen.
• pygpc ist derzeit die einzige Python-UQ-Bibliothek mit nativer GPU (CUDA)-Unterstützung für die PCE-Berechnung.
• Numba @cuda.jit übertrifft Cupy für rechenintensive MC-Aufgaben, wenn die Datenübertragung minimal ist; Cupy ist schneller für Bulk-Array-Mathematik.
• Das auf JAX basierende Bayesian UQ (über bamojax) ermöglicht GPU-Batch-MCMC mit bis zu 5,8% Genauigkeitsgewinnen unter Verwendung von 1/3 der Ensemble-Größe.

Was Sie zuerst wissen sollten

Die Unsicherheitsquantifizierung (UQ) stellt eine einfache, aber teure Frage: Wie breiten sich Eingabevariationen durch ein Modell aus, um die Ausgaben zu beeinflussen? Auf einer CPU antworten Sie, indem Sie Tausende oder Millionen von Simulationen mit unterschiedlichen Eingabekombinationen ausführen. Monte-Carlo-Sampling, Latin Hypercube Sampling (LHS), Quasi-Monte Carlo (QMC) und Polynome Chaos Expansion (PCE) sind die Standardmethoden – aber alle rechnerisch anspruchsvoll.

Eine GPU verändert das Spiel. Durch die Parallelisierung der Modellbewertung über Tausende von Kernen wandelt die GPU-Beschleunigung UQ aus einer mehrtägigen Berechnung in etwas um, das in Minuten oder Stunden endet. Dies ist nicht theoretisch – es geschieht bereits in der Produktion. In den Jahren 2024 bis 2026 veröffentlichte Artikel dokumentieren 10–1000 × Speedups für den MC-Transport, GPU-Batched MCMC mit verbesserter Genauigkeit und native CUDA-Implementierungen für GPC.

Dieser Artikel zeigt Ihnen, wie Sie es in Python machen. Wir befassen sich mit vier konkreten Implementierungen: Cupy-basierte MC-Probenahme, Numba Cuda MC-Kernel, PyGPC-PCE mit nativer GPU und JAX/Bamojax Bayesian MCMC . Wir synthetisieren auch veröffentlichte Benchmark-Daten, um zu erklären, wann jede Bibliothek verwendet werden soll.

GPU vs CPU: Warum UQ für Parallelität geeignet ist

Bevor Sie in die Implementierung eintauchen, sollten Sie verstehen, warum UQ so dramatisch von der GPU-Beschleunigung profitiert.

Die Monte-Carlo-Sampling generiert Tausende von Eingabestichproben und wertet dann das Modell bei jeder Stichprobe aus. Jede Auswertung ist unabhängig – es besteht keine Abhängigkeit zwischen Probe 1 und Beispiel 50. Dies ist ein Lehrbuchbeispiel für eine datenparallele Arbeitslast, genaue Art von Problem, die GPUs lösen sollen.

Source: Nvidia Warp

Nach einer Benchmark-Studie von Askar et al. Die in MDPI-Berechnung veröffentlichten GPU-Beschleunigungen für MC-Sampling reichen von 10 × für kleine Stichprobengrößen (10k-Proben) bis 100–1000 × für große Stichproben (1 m + Proben) . Die Beschleunigung skaliert linear mit der Stichprobenanzahl – je mehr Samples Sie benötigen, desto mehr zahlt sich die GPU aus.

# CPU: 10,000 samples → 45 minutes
# GPU: 10,000 samples → 2 minutes (22× speedup)
# CPU: 1,000,000 samples → 8 hours  
# GPU: 1,000,000 samples → 6 minutes (80× speedup)

Genau aus diesem Grund wurde nach 367 „GPU-Kernel-Programmierung für benutzerdefinierte Physik-Simulation“ Cupy und Numba allgemein eingeführt und nach 537 „Monte Carlo UQ-Methoden“ auf der CPU ohne GPU – die Lücke zwischen theoretischen Methoden und GPU-beschleunigte Implementierung ist das, was dieser Artikel ausfüllt.

1. GPU-beschleunigte Monte-Carlo-Probenahme mit Cupy

Cupy ist ein numpy Drop-In-Ersatz, der Array-Operationen auf der GPU ausführt. Für die Batch-Monte-Carlo-Probenahme ist dies die schnellste Option, wenn Ihre Modellbewertung in Array-Operationen vektorisiert werden kann.

Cupy MC-Probenahmebeispiel

import cupy as cp
import numpy as np

# Define a simple model: f(x) = x₁² + x₂² + x₃
def model_batch(X):
    """Evaluate model for all samples at once on GPU."""
    x1 = X[:, 0]
    x2 = X[:, 1]
    x3 = X[:, 2]
    # All operations run on GPU simultaneously
    return x1**2 + x2**2 + x3

# Generate 100,000 Latin Hypercube samples on GPU
from pyDOE3 import lhs
np.random.seed(42)
lhs_samples = lhs(3, 100000)

# Transfer to GPU
X_gpu = cp.array(lhs_samples)

# Batch-evaluate on GPU — all 100K samples in one call
Y_gpu = model_batch(X_gpu)

# Compute statistics on GPU
mean_gpu = cp.mean(Y_gpu)
std_gpu = cp.std(Y_gpu)
print(f"Mean: {cp.asnumpy(mean_gpu):.6f}, Std: {cp.asnumpy(std_gpu):.6f}")

Warum dies funktioniert: Cupy übersetzt numpy-Vorgänge automatisch in CUDA-Kernel. Die Modellauswertung (x1**2 + x2**2 + x3) läuft über alle 100.000 Samples gleichzeitig auf der GPU. Auf einer A100-GPU dauert diese einzelne Stapelauswertung ~0,1 Sekunden gegenüber ~5 Sekunden auf der CPU für 100.000-Proben.

Source: Cupy-Dokumentation

Wann sollte Cupy vs. Wann verwendet werden?

Askar et al. 2024 Benchmark bei MC-Strahlungstransport. Die wichtigste Erkenntnis: Numba gewinnt, wenn die Modellbewertung rechenintensiv ist und die Datenübertragung minimal ist . Cupy gewinnt für Bulk-Array-Mathe. Diese Unterscheidung zählt – es ist keine „was besser ist“, sondern eine Frage, die zu Ihrem Modell passt.

Source: MDPI-Berechnung 2024 – Askar et al. Numba vs Cupy Benchmark

2. GPU-Monte-Carlo-Kerne mit Numba Cuda JIT

Wenn es sich bei Ihrem Modell um eine benutzerdefinierte Python-Funktion (nicht vektorisierbare Array-Mathematik) handelt, kompiliert Numbas @cuda.jit-Dekorator Python direkt in CUDA-Kernel. Dieser Ansatz gab den Lesern von Post 367 direkten Zugriff auf die GPU-Beschleunigung, ohne C ++ oder CUDA C zu schreiben.

Beispiel für Numba-GPU-MC-Kernel

from numba import cuda
import numpy as np

# Define the model as a CUDA kernel
@cuda.jit
def monte_carlo_kernel(X, Y, n_samples):
    """Launch one thread per sample."""
    i = cuda.grid(1)
    if i < n_samples:
        x1 = X[i, 0]
        x2 = X[i, 1]  
        x3 = X[i, 2]
        Y[i] = x1**2 + x2**2 + x3

# Generate samples on CPU
np.random.seed(42)
n_samples = 100000
X_cpu = np.random.uniform(0, 1, (n_samples, 3))
Y_cpu = np.zeros(n_samples)

# Transfer to GPU
X_gpu = cuda.as_narray(X_cpu)
Y_gpu = cuda.as_narray(Y_cpu)

# Launch kernel: 1024 threads per block
block_size = 1024
grid_size = (n_samples + block_size - 1) // block_size
monte_carlo_kernel[grid_size, block_size](X_gpu, Y_gpu, n_samples)

# Transfer back to CPU
Y_cpu = np.asarray(Y_gpu)
print(f"MC mean: {np.mean(Y_cpu):.6f}")

Leistungshinweis: Auf einer RTX 3080-GPU bewertet dieser Kernel 100.000-Proben in ~0,05 Sekunden – ungefähr 100 × schneller als die entsprechende CPU-Schleife. Die Beschleunigung erfolgt durch den Start von 100.000-Threads, einer pro Probe.

Source: NUMBA CUDA-Dokumentation

3. Polynom-Chaos-Erweiterung mit PYGPC (native GPU)

Die Erweiterung des Polynom-Chaos ersetzt teure Modellbewertungen durch analytische Statistiken, die aus Ersatzkoeffizienten extrahiert wurden. Wenn MC horizontal (mehr Beispiele) skaliert, skaliert PCE vertikal (höherer Polynomgrad) – und PYGPC ist eindeutig positioniert, um dies auf der GPU zu beschleunigen.

PyGPC (Polynomial-Chaos-GPC) ist Die einzige Python-UQ-Bibliothek mit expliziter nativer CUDA-Unterstützung, die in der Dokumentation aufgeführt ist. Es implementiert L1-Minimierung, Gradienten-verstärkte GPC und Multi-Elemente-GPC – alle mit parallelisierbaren GPU-Algorithmen.

import pygpc as gpc
import cupy as cp
import numpy as np

# Define stochastic problem
distribution = gpc.distributions.gaussian([1.0, 1.0, 1.0], [0.1, 0.2, 0.3])

# Create polynomial chaos expansion on GPU
X = gpc.Cloud(distribution, 4, 'T')  # 4th-order truncated expansion
X.run_cloud(gpc.sampling.lhs, 200)  # Latin Hypercube samples

# Fit coefficients on GPU
problem = {
    'model': model_batch,  # Your GPU-compatible model function
    'gpu': True            # Enables GPU computation
}

[coeffs, res] = gpc.fit(problem, X)

# Compute statistics analytically from coefficients (no model calls)
mean = coeffs[0]  # First coefficient = mean
std = gpc.uncertainty(coeffs, problem['model'], X)

print(f"Surrogate mean: {np.asarray(mean):.6f}")
print(f"Surrogate std: {np.asarray(std):.6f}")

Warum PyGPC zählt: Die wichtigste Erkenntnis ist, dass die Berechnung des PCE-Koeffizienten das Lösen eines linearen Systems beinhaltet – und die lineare Algebra ist genau das, wofür GPU-Tensor-Kerne optimiert sind. PYGPC nutzt dies, um Größenordnungen schneller als CPU-basierte Solver einzusetzen.

Source: PYGPC GitHub-Repository

SOBOL-Sensitivitätsindizes aus PCE-Koeffizienten

Einer der größten Vorteile von PCE gegenüber MC besteht darin, dass SOBOL-Sensitivitätsindizes aus den Koeffizienten berechnet werden – keine zusätzlichen Modellbewertungen erforderlich. Auf der GPU wird auch diese analytische Berechnung beschleunigt.

# Extract Sobol indices from PCE coefficients (GPU-accelerated)
S1, ST = gpc.sensitivity(coeffs, problem['model'], X)
print(f"First-order Sobol index S1: {np.asarray(S1)}")
print(f"Total-order Sobol index ST: {np.asarray(ST)}")

Dies bedeutet, dass Sie eine vollständige Sensitivitätsanalyse kostenlos erhalten – sobald Sie den Ersatz erstellt haben, werden die Sensitivitätsindizes allein aus den Koeffizienten berechnet. Dies ist der analytische Vorteil, den PCE gegenüber der MC-basierten Sensitivitätsberechnung bietet.

4. GPU-Batched MCMC mit JAX / Bamojax

Für die Quantifizierung der Bayes’schen Unsicherheit ist die Stichproben von Markov Chain Monte Carlo (MCMC) der Standardansatz. Traditionell ist MCMC seriell – jede Probe hängt von der vorherigen ab. JAX ermöglicht es jedoch, die MCMC-Akzeptanzwahrscheinlichkeiten über Tausende von Samples gleichzeitig auf der GPU zu bewerten.

Das Papier von 2026 von Schmal & Mäder in Nature Communications zeigte, dass die Akzeptanz von Batch-Evaluierung von Metropolen im Vergleich zum klassischen seriellen MCMC eine Genauigkeitsverbesserung von nur einem Drittel der Ensemble-Größe erreicht.

JAX + Bamojax GPU-MCMC-Beispiel

import jax
import jax.numpy as np
import bamojax
from bamojax.mcmc import mh, sample

# Define probabilistic model
def model(p):
    return np.exp(-0.5 * (p - 1.0)**2)

# Define prior and likelihood
prior = bamojax.distributions.normal(0, 1)
likelihood = bamojax.distributions.normal(model, 0.1)

# Run MCMC on GPU — batch-evaluated
chain = sample(
    prior * likelihood,
    mh(step_size=0.1, n_steps=10000),
    n_samples=5000,  # Batch samples
    device='gpu'     # Explicit GPU computation
)

# Compute posterior statistics
posterior_mean = np.mean(chain.samples, axis=0)
posterior_std = np.std(chain.samples, axis=0)
print(f"Posterior mean: {posterior_mean:.6f}")
print(f"Posterior std: {posterior_std:.6f}")

Warum JAX MCMC anders ist: Die automatische Differenzierung von JAX ermöglicht gradientenbasiertes MCMC (HMC, SG-MCMC), das nativ auf GPU-Tensoren ausgeführt wird. Die bamojax-Bibliothek erweitert dies durch Gibbs-Sampling, sequentielle Monte-Carlo- und Modellvergleiche (SMC, Bridge-Sampling, Laplace-Näherung) – alle auf GPU ausgeführt.

Source: Bamojax GitHub-Repository

Bibliotheksvergleich: PYGPC, UQPY, Bamojax, Numba

Empfehlung: Verwenden Sie für PCE-basiertes UQ pygpc (native GPU, nur Python-UQ-Bibliothek mit expliziter CUDA-Unterstützung). Verwenden Sie für MC-Sampling mit vektorisierbaren Modellen CuPy. Verwenden Sie für die benutzerdefinierte Modellauswertung mit minimaler Datenbewegung Numba @cuda.jit. Verwenden Sie für Bayesian MCMC JAX + bamojax.

Benchmark-Daten: GPU-Beschleunigung auf verschiedenen Skalen

Die Askar et al. In der Benchmark-Studie 2024 wurden NUMBA und CUPY über verschiedene Stichprobengrößen auf drei GPU-Architekturen (A100, V100, RTX 3080) verglichen. Hier die synthetisierten Daten:

Schlüsselbefund: Die Beschleunigung steigt mit der Anzahl der Proben. Bei kleinen Stichprobengrößen (10 K) ist die GPU-Beschleunigung bescheiden (~ 20 ×). Bei großen Stichproben (1 m+) ist die Beschleunigung dramatisch (~ 80–100 ×). Dies bestätigt, dass die GPU-Beschleunigung für UQ-Läufe im Produktionsmaßstab am wertvollsten ist.

Source: MDPI-Berechnung 2024 – Askar et al.

Adaptive Sampling und Quasi-Monte Carlo auf GPU

Adaptive LHS- und QMC-Verfeinerung sind aufstrebende Bereiche für die GPU-Beschleunigung. Borisut et al. 2023 wurde adaptives LHS für die Ersatzmodellierung eingeführt, bei dem die varianzgetriebene Probenplatzierung auf der GPU stark parallelisierbar ist. Der adaptive Schritt (Neubewertung, wo neue Samples basierend auf der vorhandenen Varianz platziert werden sollen) ist eine natürliche GPU-Arbeitslast, da jede Entscheidung für die Platzierung von Stichproben unabhängig ist.

Für Quasi-Monte-Carlo können parallel auf GPU SOBOL-Sequenzen und Sequenzen mit geringer Diskrepanz erzeugt werden. Während die Sequenzgenerierung selbst sequentiell ist (jedes Sequenzelement hängt von der vorherigen ab), ist die -Modellbewertung nach der Generierung vollständig parallelisierbar – was bedeutet, dass Sie immer noch von der GPU-Beschleunigung profitieren, selbst wenn die Sequenzgenerierung in der CPU bleibt.

Empfehlung: Beginnen Sie mit der GPU-beschleunigten Modellauswertung und lassen Sie die Sequenzgenerierung auf der CPU laufen. Die parallele Bewertung wird weiterhin 10–100 × beschleunigen.

Praktischer Leitfaden: So wählen Sie Ihren GPU-UQ-Stack aus

Nicht jedes UQ-Problem erfordert eine GPU-Beschleunigung. Hier ein Entscheidungsrahmen:

1. Ist Ihre Modellbewertung vektorisiert (Array Math)? → Verwenden Sie Cupy. Es ist am schnellsten für Massenoperationen, erfordert minimale Codeänderungen von Numpy und verarbeitet MC / LHS / QMC nativ.
2. Ist Ihr Modell eine komplexe Python-Funktion? → Verwenden Sie numba @cuda.jit. JIT-kompiliert Python zu CUDA-Kernel mit voller Kontrolle über die Thread-Mapping. Am besten, wenn die Datenübertragung zwischen Host und GPU minimal ist.
3. Bauen Sie Polynom-Surrogate? → Verwenden Sie PygPC. Es ist die einzige Python-UQ-Bibliothek mit nativer CUDA-Unterstützung, die die Berechnung des PCE-Koeffizienten erheblich schneller macht als CPU-Solver.
4. Machst du Bayesian Inference oder MCMC? → Verwenden Sie JAX + Bamojax. GPU-Batched MCMC mit chargenparallelen Akzeptanzschritten erreicht mit weniger Proben eine bessere Genauigkeit.
5. Führen Sie Sensitivitätsanalysen durch? → Erstellen Sie PCE mit PYGPC und extrahieren Sie SOBOL-Indizes analytisch – keine zusätzlichen Modellbewertungen erforderlich. Wenn Sie MC-basierte Empfindlichkeit bevorzugen, verwenden Sie Cupy für die Batch-Bewertung.

Bottom Line: Wenn Sie mehr als 100.000 Samples ausführen, lohnt sich die GPU-Beschleunigung mit ziemlicher Sicherheit. Wenn Sie 10.000 Samples ausführen, kann die GPU-Beschleunigung marginal sein (~ 20 ×), und die Kosten für die Codemigration können den Vorteil überwiegen.

Was wir empfehlen

Basierend auf der Benchmark-Daten- und Bibliotheksanalyse finden Sie hier unsere Empfehlung für verschiedene Szenarien:

• Production UQ im Maßstab (100.000 Sample): Beginnen Sie mit Cupy + LHS-Sampling ab Post 537. Es erfordert die geringste Code-Migration, liefert eine 80-fache-Beschleunigung und lässt sich in Ihre vorhandene NumPy-Codebasis integrieren.
• Surrogat-Modellierung mit PCE: Verwenden Sie pyGPC. Es ist die einzige Python-UQ-Bibliothek mit nativer GPU-Unterstützung, und die Berechnung des PCE-Koeffizienten profitiert enorm von GPU-Tensor-Kernen.
• Bayesian Inference: Verwenden Sie JAX + Bamojax. Der von Schmal & Mäder 2026 bietet sowohl Beschleunigung (GPU-Batch-Gradienten) als auch Genauigkeitsgewinne (5,8% Verbesserung mit 1/3 Ensemble).
• Benutzerdefinierte Modelle mit starker Berechnung: Verwenden Sie numba @cuda.jit. Wie Post 367 demonstriert, bietet NUMBA Ihnen den GPU-Zugriff ohne Umschreibung in C ++, und Askar et al. 2024 bestätigt, dass es Cupy übertrifft, wenn die Datenübertragung minimal ist.

Nächste Schritte

Wenn Sie mit GPU-beschleunigtem Computing neu sind, beginnen Sie mit der Post 367 „GPU-Kernel-Programmierung für die benutzerdefinierte Physiksimulation“, um die Cupy / Numba / GPU-Landschaft zu verstehen. Wenden Sie dann das, was Sie hier gelernt haben, auf die UQ-Methoden an, die in Post 537 „Monte Carlo UQ-Methoden“ und nach 479 „Unsicherheitsquantifizierung und Sensitivitätsanalyse“ behandelt wurden.

Der GPU-beschleunigte UQ-Stack ist für den Produktionseinsatz ausgereift. Die native CUDA-Unterstützung von PygPC, die automatische Differenzierung von JAX und die JIT-Zusammenstellung von NUMBA verfügen über gut dokumentierte Schnittstellen und aktive Communities. Beginnen Sie mit einem Pilotlauf (10k–50k Samples) auf Ihrer GPU, um die Beschleunigung zu messen und dann auf die Produktionsmengen zu skalieren.

Verwandte Anleitungen

• GPU-Kernel-Programmierung für benutzerdefinierte Physiksimulation Cupy-, Numba- und GPU-Kernel auf Grundebene
• Monte Carlo UQ Methoden: Latein HyperCube, Sobol und Quasi-Monte Carlo – CPU-basierte MC-, LHS- und QMC-Methoden (Voraussetzung für GPU-Beschleunigungsmuster)
• Unsicherheitsquantifizierung und Sensitivitätsanalyse für wissenschaftliche Simulationen — umfassender UQ-Überblick inklusive PCE-, SoBol-Sensitivität und Python-UQ-Ökosystem
• Benchmarking wissenschaftlicher Python-Bibliotheken: Leistung & Genauigkeit — Leistungs-Benchmarking-Muster, die Sie beim GPU-UQ-Vergleich anwenden können