POSSIS: 3D Radiative Transfer in Transients

• POSSIS is a Monte Carlo radiative-transfer code for explosive transients that computes angle-dependent spectra, light curves, and polarization with time- and wavelength-sensitive opacities.
• It incorporates detailed ejecta properties and local physics to simulate anisotropic emission, making it vital for interpreting events like GW170817 and AT2017gfo.
• POSSIS supports advanced inference through neural-network surrogates and simulation-based emulators, enhancing parameter estimation and survey design.

POSSIS, short for POlarization Spectral Synthesis In Supernovae, is a time-dependent, three-dimensional Monte Carlo radiative-transfer code developed to predict viewing-angle dependent spectra, light curves, and polarization for multi-dimensional models of supernovae and kilonovae. In its original formulation, it propagates photon packets through homologously expanding ejecta on a 3D Cartesian grid using wavelength- and time-dependent opacities supplied as inputs; later implementations introduced nuclear heating rates, thermalization efficiencies, and wavelength-dependent opacities that depend on local ejecta properties and time. Within kilonova studies, POSSIS is used both as a direct forward model and as the physical basis for emulators, Bayesian inference pipelines, and observing-strategy studies (Bulla, 2019, Bulla, 2022).

1. Origin and domain of application

POSSIS was introduced by Bulla as a general Monte Carlo transport framework for multi-dimensional explosive transients, with explicit support for supernovae and kilonovae. Its defining capability is the production of angle-dependent observables from anisotropic ejecta, rather than angle-averaged photometry alone. This feature made it immediately relevant to the interpretation of GW170817 / AT2017gfo, where geometry, composition, and line of sight all affect the inferred ejecta properties (Bulla, 2019).

The code’s scope later expanded in two directions. First, kilonova applications became substantially more detailed, incorporating local-property-dependent heating, thermalization, and opacity prescriptions, and enabling systematic studies of how simplified microphysics biases spectra and light curves (Bulla, 2022). Second, the same transport machinery was adapted outside the kilonova context. In modeling tidal disruption events in the collision-induced outflow scenario, POSSIS was used in a grey, Thomson-scattering configuration to compute continuum polarization as a function of viewing angle and optical depth structure, showing that its packet-based formalism is not restricted to merger ejecta (Charalampopoulos et al., 2022).

This breadth of use suggests that POSSIS is best understood not as a single kilonova template, but as a radiative-transfer platform whose astrophysical specialization is determined by the adopted density, composition, opacity, and source prescriptions.

2. Radiative-transfer formulation

In the original implementation, POSSIS operates on a 3D Cartesian grid specified at a reference time $t_0$. Each cell carries velocity, density, temperature, and, in kilonova applications, electron fraction. Homologous expansion is assumed, with

$\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$

and

$T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$

Photon packets are then propagated in the lab frame while interactions are treated in the comoving frame. Initial packet frequencies are sampled from the thermal emissivity

$S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$

and interaction distances are sampled from optical depth, with the homogeneous-segment form

$s=-\frac{\ln\xi}{\kappa\rho}.$

Supported opacity channels include bound-bound, electron-scattering, bound-free, and free-free contributions; line transport can be handled either in a Sobolev mode or in an absorption-mode pseudo-continuum treatment. Re-emission uses a two-level atom scheme with redistribution probability $\epsilon=0.9$ (Bulla, 2019).

Polarization is tracked through Stokes parameters. Packets are initialized unpolarized and updated at electron-scattering events using Thomson-scattering Mueller matrices. The emergent linear polarization is reported through

$P=\frac{\sqrt{Q^2+U^2}}{I}.$

For flux extraction, POSSIS uses both direct counting and event-based estimators; the latter spawns virtual packets toward predefined observer directions and weights them by the probability of escape, substantially reducing Monte Carlo noise (Bulla, 2019).

A major later development was the replacement of spatially uniform or grey prescriptions by ejecta-local, time-dependent physics. In the improved version described in 2022, nuclear heating rates, thermalization efficiencies, and wavelength-dependent opacities depend on local ejecta properties and time. That study compared a $\texttt{FIDUCIAL}$ setup with $\texttt{SIMPLE-HEAT}$, $\texttt{SIMPLE-THERM}$, $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$0, and $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$1 variants, and found deviations of several $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$2–$\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$3 magnitudes, generally largest for $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$4 and $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$5 (Bulla, 2022).

3. Ejecta parameterizations and kilonova model families

The kilonova configuration introduced in the original POSSIS study is axisymmetric and two-component. A lanthanide-rich equatorial component spans $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$6 around the merger plane, while a lanthanide-free polar component occupies higher latitudes. In that setup, the ejecta density at $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$7 follows $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$8 with $\rho(t)=\rho_0\left(\frac{t}{t_0}\right)^{-3},$9, the temperature is initialized at $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$0, and the total heating injected into thermal bands follows Korobkin et al. with $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$1. A model with total ejecta mass $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$2 and half-opening angle $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$3 was found to provide a good match to GW170817 / AT2017gfo for orientations near the polar axis (Bulla, 2019).

Later work broadened the parameterizations. In the neural-network inference study of Pang et al., POSSIS was used in an updated “possis2” configuration with a two-component SSCr geometry: a slower post-merger wind and a faster dynamical ejecta component, including reprocessing between them. The density law is homologous with $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$4, $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$5, $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$6, and $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$7. The angular dependence of the dynamical ejecta was tested both as $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$8 and in a logistic-transition “Kawaguchi model” parameterized by $T(t)=T_0\left(\frac{t}{t_0}\right)^{-\alpha}.$9. Their 5D grids spanned $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$0, $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$1, $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$2, $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$3, and $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$4, and angular dependence improved the fit to the AT2017gfo light curves by up to $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$5 in weighted Mean Square Error (Almualla et al., 2021).

Within NMMA, the widely used Bu2019 POSSIS grid adopts a different abstraction: a two-component, uniform-composition ejecta with a lanthanide-rich equatorial region and a lanthanide-poor polar region. The boundary is set by a half-opening angle measured from the orbital plane. In the comparison to SuperNu, this Bu2019 grid is notable for not varying ejecta velocities, so posteriors are displayed using the average velocities used to expand POSSIS materials rather than velocity posteriors themselves (King et al., 22 May 2025).

POSSIS has also been specialized to black hole–neutron star mergers. In that application, the wind ejecta are spherically symmetric with $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$6, the dynamical ejecta are concentrated toward the merger plane with $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$7, the dynamical component uses $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$8, and the wind uses a flat distribution in $S(\nu)=B(\nu,T)\kappa_{\rm tot}(\nu),$9 between $s=-\frac{\ln\xi}{\kappa\rho}.$0 and $s=-\frac{\ln\xi}{\kappa\rho}.$1. Simulations were run with $s=-\frac{\ln\xi}{\kappa\rho}.$2 Monte Carlo quanta on a $s=-\frac{\ln\xi}{\kappa\rho}.$3 Cartesian grid, for 11 viewing angles and 100 logarithmically spaced epochs from 0.1 to 30 days, although only $s=-\frac{\ln\xi}{\kappa\rho}.$4–30 day light curves were analyzed because the adopted opacities are likely underestimated at earlier times (Mathias et al., 2023).

4. Observables, emulators, and inference pipelines

POSSIS natively produces time-dependent spectra and broad-band light curves, and in the original kilonova calculations synthetic observables were computed for $s=-\frac{\ln\xi}{\kappa\rho}.$5 viewing angles equally spaced in $s=-\frac{\ln\xi}{\kappa\rho}.$6 from pole to equator (Bulla, 2019). In later kilonova inference work, output bands commonly include $s=-\frac{\ln\xi}{\kappa\rho}.$7, and extended analyses add filter-specific products for PS1, 2MASS, Rubin, and JWST/MIRI (Almualla et al., 2021, King et al., 22 May 2025).

Because direct Monte Carlo evaluation is expensive, a substantial computational ecosystem has grown around POSSIS. Pang et al. built a neural-network surrogate by applying SVD to the light-curve matrices, truncating to 100 basis vectors, and training five-layer fully connected ReLU networks, one per photometric band and one for bolometric luminosity. Training all ten networks took $s=-\frac{\ln\xi}{\kappa\rho}.$8 minutes; the surrogate loaded in $s=-\frac{\ln\xi}{\kappa\rho}.$9 across all bands, compared with $\epsilon=0.9$0 minutes for a Gaussian Process Regression interpolator, and overall parameter inference took $\epsilon=0.9$1 minutes excluding loading (Almualla et al., 2021).

NMMA uses such emulators operationally. In the Bu2019lm surrogate described for low-latency gravitational-wave alerts, NMMA starts from a grid of POSSIS simulations across ejecta parameters and uses a neural network to interpolate to arbitrary inputs. The resulting products are angle-dependent absolute AB magnitudes in $\epsilon=0.9$2 and $\epsilon=0.9$3, together with the derived alert quantity

$\epsilon=0.9$4

designed for integration into the IGWN low-latency infrastructure (Toivonen et al., 2024).

A more recent route is spectral emulation. The JAX-based package fiesta trains a POSSIS surrogate on $\epsilon=0.9$5 spectral-flux snapshots generated with $\epsilon=0.9$6 photon packets per run. Its default MLP FluxModel learns PCA coefficients of $\epsilon=0.9$7 over 0.2–26 days and $\epsilon=0.9$8–$\epsilon=0.9$9, enabling bandpass integration through arbitrary filters at inference time. On an NVIDIA H100 GPU, the MLP FluxModel for POSSIS trains in $P=\frac{\sqrt{Q^2+U^2}}{I}.$0 minutes, and typical absolute-magnitude mismatches are within $P=\frac{\sqrt{Q^2+U^2}}{I}.$1 mag during the bright early phases (Koehn et al., 18 Jul 2025).

A parallel development replaces explicit likelihoods altogether. A simulation-based inference framework trained a Gaussian-process emulator on $P=\frac{\sqrt{Q^2+U^2}}{I}.$2 POSSIS simulations and then learned the posterior directly from forward simulations that include emulator uncertainty. Once trained, that SBI framework generated $P=\frac{\sqrt{Q^2+U^2}}{I}.$3 posterior samples in seconds, whereas the corresponding MCMC analysis took $P=\frac{\sqrt{Q^2+U^2}}{I}.$4 hours on the same hardware (Brown et al., 13 May 2026).

5. Scientific applications

POSSIS has been central to the modeling of AT2017gfo. The original Bulla model identified a two-component ejecta with $P=\frac{\sqrt{Q^2+U^2}}{I}.$5 and $P=\frac{\sqrt{Q^2+U^2}}{I}.$6 viewed close to face-on as a good description of the event (Bulla, 2019). Later comparative inference found that both the Bu2019 POSSIS models and a new SuperNu-based grid provide similar fits to the photometric observations, with broad agreement in wind velocity, inclination angle, and dynamical ejecta mass. The main difference is the inferred wind mass: the Bu2019 POSSIS grid returns $P=\frac{\sqrt{Q^2+U^2}}{I}.$7, whereas the preferred SuperNu TP2 morphology yields $P=\frac{\sqrt{Q^2+U^2}}{I}.$8–$P=\frac{\sqrt{Q^2+U^2}}{I}.$9; this motivated a slight preference for SuperNu in that analysis (King et al., 22 May 2025).

The code has also been used for survey and follow-up design. In optimizing optical kilonova observations with 2-m class telescopes, synthetic POSSIS light curves were combined with photon-counting noise to compare four observing sequences with the same total exposure time of 8 hours and time windows of 0.5, 1, 2, and 4 hours. That study suggested avoiding the $\texttt{FIDUCIAL}$0 filter and avoiding the use of colour curves. It also found that if the error on distance is $\texttt{FIDUCIAL}$1, the 0.5, 1, and 2-hour time-window sequences are equivalent, and recommended the 2-hour sequence because it has 1 day cadence; when the distance is unknown, the 0.5-hour sequence is preferable (Camisasca et al., 2023).

For morphology discrimination, POSSIS has served as the reference against which alternative radiative-transfer suites are judged. Rubin/LSST optical target-of-opportunity data alone were found insufficient to distinguish SuperNu morphology classes, whereas adding late-time JWST/MIRI photometry at 8, 15, and 20 days markedly improved dynamical-ejecta recovery and enabled model selection. In that context, POSSIS functions both as a competing radiative-transfer family and as a baseline for understanding which aspects of the data are actually geometry-sensitive (King et al., 22 May 2025).

In compact-binary population studies, POSSIS has been integrated with ejecta-fitting formulae and gravitational-wave posteriors. The low-latency framework based on Bu2019lm produces event-by-event $\texttt{FIDUCIAL}$2 values and predicted $\texttt{FIDUCIAL}$3 light curves, while the black hole–neutron star study of Anand et al. used POSSIS to show that a soft 2-families equation-of-state scenario should not produce a kilonova unless $\texttt{FIDUCIAL}$4, $\texttt{FIDUCIAL}$5, and $\texttt{FIDUCIAL}$6, whereas a 1-family scenario yields detectable kilonovae across a wider region of parameter space (Toivonen et al., 2024, Mathias et al., 2023).

Joint modeling has extended POSSIS beyond isolated kilonova fits. In merger-driven GRBs, NMMA combines POSSIS for the thermal kilonova component with afterglowpy for the non-thermal afterglow, fitting both simultaneously. In that sample, the median wind mass,

$\texttt{FIDUCIAL}$7

exceeded the median dynamical mass,

$\texttt{FIDUCIAL}$8

and the inferred relation

$\texttt{FIDUCIAL}$9

was reported between wind mass and beaming-corrected jet kinetic energy (Singh et al., 29 Dec 2025).

Outside kilonovae, POSSIS has been used to interpret TDE polarimetry. In the collision-induced outflow scenario, high peak fallback rates yielded polarization levels below one percent for every viewing angle, whereas low peak fallback rates produced maxima of $\texttt{SIMPLE-HEAT}$0, with the location of $\texttt{SIMPLE-HEAT}$1 shifting from equatorial to intermediate viewing angles as the optical-depth structure changed (Charalampopoulos et al., 2022).

6. Limitations, systematics, and methodological debates

A persistent feature of POSSIS is that transport is performed using prescribed opacities rather than a fully self-consistent microphysical solution. In the original code description, this design choice was presented explicitly as a means of retaining fast computations while modeling multi-dimensional geometry, anisotropic opacities, and polarization (Bulla, 2019). The consequence is that inference quality depends directly on the adopted heating, thermalization, and opacity prescriptions.

This dependence has become a major subject of scrutiny. The improved 2022 POSSIS study showed that replacing local, time-dependent treatments by uniform heating, constant thermalization, or grey opacities can shift synthetic observables by several $\texttt{SIMPLE-HEAT}$2–$\texttt{SIMPLE-HEAT}$3 magnitudes, with the largest discrepancies in the $\texttt{SIMPLE-HEAT}$4 and $\texttt{SIMPLE-HEAT}$5 models (Bulla, 2022). In the neural-network AT2017gfo study, LTE was explicitly noted as likely to fail at late times when the ejecta become optically thin, motivating early-time weighting in the likelihood (Almualla et al., 2021). In the BH–NS calculations, epochs earlier than 0.5 days were excluded because the adopted opacities are likely underestimated there (Mathias et al., 2023).

Model geometry is another systematic. The Bu2019 POSSIS grid used in several NMMA analyses does not vary ejecta velocities, and its geometry is encoded through a two-component uniform-composition model with a half-opening angle, whereas the SuperNu comparison study explored multiple axisymmetric morphologies with 54 angular bins. That difference was identified as one reason why the two radiative-transfer families can return different ejecta parameters even when they provide comparably good photometric fits (King et al., 22 May 2025).

The growth of emulation has introduced a further layer of approximation. The SBI analysis based on a Gaussian-process POSSIS emulator found that MCMC with a diagonal Gaussian likelihood can suffer from systematic bias because the likelihood fails to capture the non-Gaussian, correlated structure of emulator uncertainty. In that analysis, some MCMC posteriors piled up at prior boundaries, while SBI learned the posterior directly from forward simulations that included the full predictive distribution (Brown et al., 13 May 2026). A plausible implication is that current high-throughput use of POSSIS increasingly depends not only on the fidelity of the radiative-transfer calculation itself, but also on how emulator uncertainty is represented.

Taken together, these issues define the present status of POSSIS. It is a widely used forward model for anisotropic transient emission, especially kilonovae, but its scientific output is inseparable from assumptions about ejecta morphology, opacity data, thermalization, LTE validity, and statistical treatment of emulator error.