#  WHOLE DEVICE MODELING OF THE FUZE SHEARED-FLOW-STABILIZED Z
PINCH

I. A. M. DATTA

Zap Energy

Seattle, USA

Email<EMAIL_ADDRESS>
E. T. MEIER

Zap Energy

Seattle, USA

U. SHUMLAK

Zap Energy / University of Washington

Seattle, USA

## Abstract

The FuZE sheared-flow-stabilized Z pinch at Zap Energy is simulated using
whole-device modeling employing an axisymmetric resistive magnetohydrodynamic
formulation implemented within the discontinuous Galerkin WARPXM framework.
Simulations show formation of Z pinches with densities of approximately 1022
m-3 and total DD fusion neutron rate of 107 per $\mu$s for approximately 2
$\mu$s. Simulation-derived synthetic diagnostics show peak currents and
voltages within 10% and total yield within approximately 30% of experiment for
similar plasma mass. The simulations provide insight into the plasma dynamics
in the experiment and enable a predictive capability for exploring design
changes on devices built at Zap Energy.

## 1\. Introduction

Whole device modeling of the FuZE shear-flow-stabilized (SFS) Z pinch [1, 2]
at Zap Energy is illuminating key physics of the device and guiding progress
towards a breakeven fusion concept with a high-gain power-producing power
plant. Two-dimensional axisymmetric resistive magnetohydrodynamic (MHD)
simulations using the discontinuous Galerkin finite element WARPXM framework
[3] are compared to experimental measurements from the FuZE device using a
suite of synthetic diagnostics. The modeling unveils the plasma dynamics
inside FuZE and enables in-depth investigation of various areas including
plasma compression adiabaticity and electrode erosion. Validated modeling
provides a predictive tool for next-step SFS Z pinch experimental design.

This paper summarizes simulation results using synthetic diagnostics to
compare with experimental analogues. Section 2 describes the MHD model used
while Sec. 3 briefly describes the numerical method. Section 4 describes the
problem setup, Sec. 5 discusses results, and Sec. 6 gives conclusions.

## 2\. Model

For FuZE plasmas, expected temperatures are on the order of 1 keV for a pinch
size of about 1 cm with enclosed currents of up to 500 kA [1, 4]. In this
regime, the ion Larmor radius and ion inertial length are on the millimeter
scale. For the MHD model to be valid, these length scales should be small
compared to the system [5, 6]. For a 1 keV pinch with 500 kA current, a radius
of 1 cm would be on the order of 10 ion Larmor radii, so we proceed with a
resistive MHD model given by

$\frac{\partial\rho}{\partial
t}+\nabla\cdot\left(\rho\mathbf{u}\right)=0\mathrm{,}$ (1)
$\frac{\partial\left(\rho\mathbf{u}\right)}{\partial
t}+\nabla\cdot\left(\rho\mathbf{u}\mathbf{u}+p\overline{\overline{I}}+\overline{\overline{\Pi}}\right)=\left(\omega_{c}\tau\right)\mathbf{J}\times\mathbf{B}\mathrm{,}$
(2) $\frac{\partial e_{t}}{\partial
t}+\nabla\cdot\left[\left(e+p\right)\mathbf{u}+\overline{\overline{\Pi}}\cdot\mathbf{u}+\mathbf{h}\right]=-\nabla\cdot\left(\mathbf{E}\times\mathbf{B}\right)-C_{bb}\frac{\rho}{A_{\mathrm{i}}}T^{4}\mathrm{,}$
(3) $\frac{\partial\mathbf{B}}{\partial t}+\nabla\times\mathbf{E}=0\mathrm{,}$
(4)
$e_{t}=e+\frac{\mathbf{B}\cdot\mathbf{B}}{2}=\frac{p}{\gamma-1}+\frac{\rho\mathbf{u}\mathbf{u}}{2}+\frac{\mathbf{B}\cdot\mathbf{B}}{2}\mathrm{,}$
(5)
$\mathbf{J}=\frac{\nabla\times\mathbf{B}}{\left(\omega_{c}\tau\right)}\mathrm{,}$
(6)
$\mathbf{E}=-\mathbf{u}\times\mathbf{B}+\frac{\left(\nu_{p}\tau\right)}{\left(\omega_{c}\tau\right)}\eta\mathbf{J}\equiv\mathbf{E}_{\mathrm{Ideal}}+\mathbf{E}_{\mathrm{Resistive}}\mathrm{,}$
(7)
$\overline{\overline{\Pi}}=-\mu\left[\nabla\mathbf{u}+\left(\nabla\mathbf{u}\right)^{T}-\frac{2}{3}\left(\nabla\cdot\mathbf{u}\right)\overline{\overline{I}}\right]\mathrm{,}$
(8) $\mathbf{h}=-k\nabla T\mathrm{.}$ (9)

Equations (1)-(9) are written in non-dimensional form with normalized
cyclotron frequency $\left(\omega_{c}\tau\right)$ and collision frequency
$\left(\nu_{p}\tau\right)$ defined as described in Ref. [7]. The equations are
solved, as described in the next section, in two-dimensional ($r$-$z$ plane)
axisymmetric form where $\frac{\partial}{\partial\theta}=0$.

A blackbody radiation model [8] is added in Eq. (3). Viscosity and thermal
conduction terms are added as well, as seen in Eqs. (2), (3), (8), and (9).
These terms are required to match the simulation dynamics with experiment, as
described further in Sec. 5.

## 3\. Numerical Method

The MHD model described in Eqs. (1)-(9) can be combined and expressed
compactly as

$\frac{\partial\mathbf{q}}{\partial
t}+\nabla\cdot\overline{\overline{F}}=\mathbf{S}\mathrm{.}$ (10)

The equations are solved using the discontinuous Galerkin method within the
WARPXM framework developed at the University of Washington and Zap Energy [3].
An overview of the method and implementation is found in Refs. [3], [5], and
[7]. The discontinuous Galerkin method is well suited for solving the
predominantly hyperbolic MHD equations on unstructured meshes using simplex
finite elements, allowing for simulation of the FuZE geometry. In this work, a
cylindrical formulation of the MHD equations described in Sec. 2 is used, the
numerical handling of which is described in Ref. [5]. Artificial dissipation
is also added to stabilize solutions [9].

A circuit model is implemented to capture the interaction between the driving
current bank and the plasma in FuZE in a self-consistent manner. The circuit
is modeled by including the plasma dynamics as a time-dependent impedance and
associated load voltage as a circuit element. The current is driven by the
circuit with prescribed capacitance, inductance, resistance, charge voltage,
and trigger time.

## 4\. Problem Setup

The simulations performed aim to match conditions for a set of experimental
operating parameters on FuZE as shown in Fig. 1. Triangular finite elements
are used with a quadratic solution basis. A slug of plasma of 0.5 mg in the
acceleration region is initialized where the valves inject gas. A volumetric
injection of plasma is also input at a rate of 3 mg / ms at the valve
location, mimicking breakdown of residual neutral gas and ongoing valve
injection.

Figure 1: Two-dimensional $r$-$z$ plane slice of FuZE used as the simulation
domain. Gas is puffed into the coaxial acceleration region where it is ionized
by an electrical power supply between the cathode (inner electrode) and anode
(outer electrode) and pushed into the assembly region to form a pinch. The
simulation assumes the puffed gas is fully ionized as an initial plasma slug
at $t=0$. Note that the radial dimension is stretched roughly by 3-to-1.

Simulated initial mass is chosen by following the trajectory of the magnetic
field data for an experimental plasma pulse (#230518072) and calculating the
corresponding plasma mass. Valve mass injection rates are also estimated, as
described in Sec. 4.1. A series RLC circuit is used to simulate the driving
capacitor bank for FuZE, with circuit values $C=222~{}\mu\mathrm{F}$,
$L=286\mathrm{~{}nH}$, and $R=1.5\mathrm{~{}m}\Omega$. At $t=0$, the capacitor
charged to 25 kV is discharged. The plasma load is calculated as a voltage
drop measured along the insulator of the device between the cathode (inner
electrode) and anode (outer electrode). The gap voltage is measured upstream
of the insulator ahead of a 70 nH region of inductance (see Fig. 1 in Ref.
[10]), so to compare with experiment, it is calculated as
$V_{\mathrm{gap}}=-\int\limits_{\text{cathode}}^{\text{anode}}E_{r}dr+\left(\mathrm{70~{}nH}\right)\dot{I}\mathrm{.}$

The circuit current, converted to an azimuthal magnetic field, imposes a
boundary condition on the plasma at the insulator, where the axial velocity is
allowed to float if into the domain and set to zero otherwise. At $r=0$,
axisymmetric boundary conditions are employed, as described in Ref. [5]. All
other boundaries are treated as free-slip conducting walls that are thermally
insulating by imposing zero-normal density and pressure gradients.

### 4.1. Mass determination from slug and valve injection models

The mass of the plasma in experiment is estimated using a combination of
magnetic probes along the outer electrode and gas valve injection formulas.
Using the magnetic probes, the motion of an assumed slug of plasma injected by
the gas puff valves is tracked to estimate its mass. Using this slug model,
the plasma mass is determined using the force [11, 12],

$F(t)=\frac{L^{{}^{\prime}}I(t)^{2}}{2}=m\ddot{z}\mathrm{.}$ (11)

Dividing Eq. (11) by $m$ and integrating twice to time, $t_{d}$, leads to an
equation for distance traveled by the plasma,

$d=\frac{L^{{}^{\prime}}}{2m}\iint\limits_{0}^{t_{d}}I^{2}dt\mathrm{.}$ (12)

In the experiment, the times taken for the plasma slug to move from the gas
puff valves to the magnetic probes along the outer electrode are known by
measuring when each probe records a jump in magnetic field. A space-time plot
of this magnetic probe data as shown in Fig. 2 indicates where and when along
the outer electrode the jumps occur. With $d$ and $t_{d}$ known for various
magnetic probe locations, $m$ is calculated using Eq. (12). The dashed orange
line in Fig. 2 corresponds to $m=0.5$ mg for the particular plasma pulse.

Figure 2: Enclosed current in experiment as measured by magnetic probes along
the outer electrode as a function of time. A few probes in the assembly region
did not register data for this pulse.

A combination of gas valve injection and a deflagration process ionizing
remaining neutral gas sustains plasma flow to the pinch [11]. To account for
that effect, a mass injection model is needed. A simple model assumes gas
flows out of the valve at the speed of sound such that

$\frac{dm}{dt}=V\frac{d\rho}{dt}=-\rho c_{s}A=-A\rho\sqrt{\frac{\gamma
p}{\rho}}\mathrm{.}$ (13)

Normalizing by non-dimensional parameters, $\tilde{\rho}=\rho/\rho_{0}$ and
$\tilde{t}=t/\tau_{v}$ where $\tau_{v}=V/c_{s0}A$ is the timescale of the
venting valve and assuming adiabatic expansion such that
$\tilde{p}=p/p_{0}=\left(\rho/\rho_{0}\right)^{\gamma}={\tilde{\rho}}^{\gamma}$
leads to
$\frac{d\tilde{\rho}}{d\tilde{t}}=-\tilde{\rho}^{\frac{\gamma+1}{2}}$. The
solution for $\tilde{\rho}\left(0\right)=1$ is
$\tilde{\rho}=\left(\frac{\gamma-1}{2}\tilde{t}+1\right)^{\frac{2}{1-\gamma}}$,
which describes the normalized density in the valve. The normalized density
entering the device is thus
$\tilde{\rho}=1-\left(\frac{\gamma-1}{2}\tilde{t}+1\right)^{\frac{2}{1-\gamma}}$.
Substitution of the dimensional values leads to an expression in terms of
initial density, $\rho_{0}$, which is recast in terms of valve pressure, $p$,
using the ideal gas law for deuterium (assuming $V=1\mathrm{~{}cm^{3}}$ and
$T=300\mathrm{~{}K}$). A reference mass calculated for a pressure of 100 PSIA
is $M_{0,\mathrm{ref}}=\rho_{0,\mathrm{ref}}V=\ 1.11\mathrm{~{}mg}$. For
different pressures, this reference mass is scaled linearly, leading to mass
as a function valve pressure and opening time as

$m=N_{v}M_{0,\mathrm{ref}}\frac{p}{100}\left[1-\left\\{1+\frac{\gamma-1}{2}\left(\frac{t-\Delta}{\tau_{v}}\right)\right\\}^{\frac{2}{1-\gamma}}\right]\mathrm{,}$
(14)

where $\Delta$ is a timing delay between the valve trigger and opening,
$N_{v}$ is the number of valves used, and $p$ is measured in PSIA. For the
experimental plasma pulse in Fig. 2, $N_{v}=4$, $p=149.1$ PSIA, and $t=0.26$
ms. For deuterium ($\gamma=7/5$, $T=300\mathrm{~{}K}$), the sound speed is 934
m/s leading to $\tau_{v}=1.070$ ms (assuming $A=1\mathrm{~{}mm^{3}}$).
Equation (14) is solved for $\Delta$ using $m=0.5$ mg, yielding $0.175$ ms.
Equation (14) is then differentiated with respect to time, yielding
$\dot{m}=3$ mg / ms. Such an injection rate may be considered a low estimate,
since deflagration is expected to make a significant contribution as an
upstream plasma source. Modeling with neutrals to capture deflagration is the
topic of ongoing research.

A simulation (#s230523033) is thus run to 30 $\mu$s matching experimental
pulse settings with 0.5 mg slug mass and 3 mg / ms injection rate. The initial
slug of plasma sits in the acceleration region where the valves are located,
as shown in Fig. 1.

## 5\. Results

Traces of current and voltage are shown in Figs. 4 and 4 with comparisons to
experimental discharge #230518072. In all experimental plots, timings are
shifted by 4.2 $\mu$s to when the current rise starts.

Figure 4 compares experimental capacitor bank current output with that of the
simulation coupling the equivalent RLC series circuit described in Sec. 4 with
the MHD plasma. The currents reach similar peaks before dropping, with
slightly lower currents later in the experimental pulse compared to
simulation.

Figure 3: Current comparison between experiment (blue) and simulation
(orange). Good agreement is found especially in earlier times.

Figure 4: $V_{\mathrm{gap}}$ comparison between experiment (blue) and
simulation (orange). Similar values and temporal profiles are seen, with jumps
early in time in simulation being an artifact of the fully-ionized MHD
assumptions.

Figure 4 shows a comparison of $V_{\mathrm{gap}}$ between experiment and
simulation. Early in time, the simulation exhibits jumps in voltage not
present in experiment, an artifact of the fully ionized MHD assumptions of the
initial slug. The incoming plasma flux is reflected by the slug, causing the
jumps in $V_{\mathrm{gap}}$. Nonetheless, the gap voltage values in experiment
and simulation are similar and have comparable temporal profiles with the drop
in voltage at around 15 $\mu$s. The voltage drop is associated with pileup of
plasma and flux at the endwall that drives a wave of magnetic flux back into
the acceleration region where it interacts with the circuit.

Figure 5 shows the evolution of plasma density along with lines of enclosed
current. The plasma slug is pushed from the acceleration region to the
assembly region through the $\mathbf{J}\times\mathbf{B}$ force imposed by the
circuit current. A true slug motion would be indicated by the lines of radial
current through the plasma as it moves through the acceleration region. At
$t=8$ $\mu$s, deviation from the 1D slug model is apparent, as axial
components of the radial contours are seen near the outer electrode. Such
behavior represents a blowby phenomenon where the imposed current does not
fully couple to the plasma [13]. The viscous ($\nu=\mu/\rho$) and heat
($\kappa=kA_{\mathrm{i}}/\rho$) diffusivities in this simulation are set to
200 m2 / s which lead to the good agreement with experimental values shown in
Figs. 4 and 4. Simulations with increased viscous and heat diffusivities show
reduced blowby and increased peak current and voltage. After the slug reaches
the assembly region and the flux rebounds after hitting the endwall, a pinch
then forms on axis, as seen in the high density on-axis regions in the 14.7
$\mu$s and 16 $\mu$s plots.

Figure 6 shows the volumetric neutron yield rate calculated using Bosch-Hale
thermonuclear fusion formulas [14] at the time of maximum total yield as well
as plasma temperature. The thermonuclear yield is a function of the density
and temperature. The combination of high temperature on axis above 3 keV
combined with on-axis densities above 1022 m-3 leads to the on-axis yield
rate, indicating the existence of a Z pinch in the simulation.

Figure 5: Plasma density at 8 $\mu$s, 12 $\mu$s, 14.7 $\mu$s (maximum total
domain-integrated yield rate time), and 16 $\mu$s, overlaid with lines of
enclosed current. The plots show current sheet canting and blowby as the
magnetic flux pushes the plasma slug into the assembly region. Plasma thermal
energy and magnetic flux accumulate at the endwall, causing a rebound behavior
that completes pinch formation. The 16-$\mu$s plot shows the pinch formed on
axis with noticeable $m=0$ activity associated with reduced sheared flow.
Figure 6: Neutron yield rate calculated from Bosch-Hale formulas for
thermonuclear fusion [14] (top) and plasma temperature (bottom) at 14.7
$\mu$s, the time of maximum total domain-integrated yield rate. Lines of
enclosed current are also shown. The plasma reaches above 3 keV near the
nosecone on axis in the assembly region. This high temperature combined with
collimated density on axis lead to the region of thermonuclear fusion seen in
the neutron yield rate plot.

Figure 7 shows the total integrated yield rate in the simulation versus
experimental data. The experimental data is of fast plastic scintillator
measurements scaled to instantaneous neutron yield rates using total yield
measurements by Lanthanum Bromide detectors. The simulation data is an
analogous synthetic diagnostic calculating the full-domain integral of the
simulated volumetric yield rate (the top plot in Fig. 6). The dashed lines are
integrations of the total yield rates over time. The difference in simulation
and experimental yield indicates an overprediction in simulation by about
300%. An adjustment to the simulation yield rate calculation is made by
removing yield where the difference between velocity at the pinch edge, the
radius $a$ where $p(a)=p_{\mathrm{max}}/4$ [5], and the axis is below 10% of
the Alfvén speed at $a$. The correction removes yield for low sheared flows
assumed to be $m=1$ unstable [15, 5]. The total yield calculated from the
correction is about 70% of the experimental yield.

In simulations without viscous and heat diffusivities (not shown), the yield
rate increases early compared to the experiment. Increased blowby in such
cases lead to earlier flux rebound from the endwall and thus earlier neutron
production. This again suggests some viscous dissipation is needed to match
experimental results. Additionally, simulations with no radiation have higher
neutron production at times after the peak, unlike in the experiment where
neutron production drops significantly after about 16 $\mu$s. The addition of
the blackbody radiation term reduces this late time yield in the simulation,
similar to experiment. The physical value for the blackbody radiation constant
used in simulations is $C_{bb}=4.559\times
10^{52}\mathrm{~{}J}^{-3}\mathrm{s}^{-1}$.

Figure 7: Total yield comparison between experiment and simulation. The
experimental yield rate (solid blue) is found using plastic scintillators
calibrated such that the time-integrated yield agrees with Lanthanum Bromide
detectors. Total volume yield rate from simulation is calculated by
integrating the volumetric yield as shown in Fig. 6 (solid orange). The solid
green line also shows the simulation volume yield rate after removing yield in
regions a pinch would be assumed unstable. Dashed lines integrate the yield
rate to determine the total yield.

Further work to simulate the plasma dynamics with more realistic viscosity and
radiation models are ongoing. This includes the consideration of Braginskii
[16] terms as well as bremsstrahlung and line-radiation models [17]. Capturing
the true nature of the endwall flux rebound is also being considered. Current
modeling does not include the possibility of outflow through the endwall which
has a spoked transparent geometry [18] or dissipation due to ion recycling and
re-ioniziation [19] that may weaken the flux rebound. However, it is unlikely
that there is full outflow [18] so, in future work, boundary conditions and
related physics (including 3D effects) will be developed to capture realistic
endwall behavior.

Figures 9 and 9 show simulation synthetic analogues of the space-time enclosed
current plot shown in Fig. 2. Figure 9 shows the enclosed current at the outer
electrode with the 0.5 mg slug trajectory superimposed. The lack of agreement
with Fig. 2 is understood by observing Fig. 5. Between 12 $\mu$s and 16
$\mu$s, a clear detachment of enclosed current is seen along the outer
electrode due to the formation of eddy currents, which reduce the total
enclosed current measured. In contrast, the experiment shows strong magnetic
fields at the outer electrode in the assembly region as seen in Fig. 2. For
the simulation, a similar plot is constructed at a smaller radius inside of
the eddy currents. Choosing 7 cm as shown in Fig. 9 indicates a magnetic
space-time structure more in line with the experiment in Fig. 2, with the 0.5
mg trajectory close to the magnetic field wave front.

Figure 8: Enclosed current in simulation at outer electrode radius. Large eddy
currents in simulation reduce these values in the assembly region.

Figure 9: Enclosed current in simulation at $r=7$ cm. At this radius, some of
the eddy currents are removed. The 0.5 mg trajectory line matches better than
in Fig. 9.

## 6\. Conclusion

Comparisons of WARPXM simulations using two-dimensional axisymmetric resistive
MHD with the FuZE experiment at Zap Energy show a remarkable level of
agreement. Comparisons using a suite of synthetic diagnostics illuminate these
agreements while suggesting areas that modeling does not fully capture.
However, these results show that the 2D axisymmetric MHD model can be used as
a predictive tool for exploring alternative geometries, electrode
configurations, driving circuit parameters, and gas fill profiles which can
accelerate progress toward breakeven fusion gain. Improvements to modeling are
possible and will improve prediction accuracy. Extension of the MHD model to
3D can capture non-axisymmetric effects leading to $m=1$ instability. The
inclusion of multi-fluid physics [20, 3, 5] can help to accurately capture
blowby, sheath effects, and instability growth rates. Multi-species continuum
kinetic modeling can allow for further accuracy of sheath modeling and pinch
dynamics at small pinch radius where the magnetic field goes to zero [21],
with hybrid approaches helping to speed up such simulations [22, 7]. Other
advancements include the addition of a plasma-neutral model [23] to capture
plasma deflagration [11] and more accurate boundary conditions that can
capture heat and particle fluxes to walls, and corresponding erosion and
impurity generation [17]. These advancements are actively being developed to
help capture the plasma dynamics of FuZE and other devices being built at Zap
Energy [4].

##  ACKNOWLEDGEMENTS

The information, data, or work presented herein was funded in part by the
Advanced Research Projects Agency – Energy (ARPA-E), U.S. Department of
Energy, under Award Nos. DE-AR-0000571, DE-AR-0001010, DE-AR-0001260 and by
the Air Force Office of Scientific Research under Grant No. FA9550-15-1-0271.
This research used resources of the National Energy Research Scientific
Computing Center (NERSC), a U.S. Department of Energy Office of Science User
Facility located at Lawrence Berkeley National Laboratory, operated under
Contract No. DE-AC02-05CH11231.

##  REFERENCES

## References

  * [1] SHUMLAK, U., Z-pinch fusion,  J. Appl. Phys. 127 (2020) 200901.
  * [2] ZHANG, Y. et al., Sustained neutron production from a sheared-flow stabilized Z pinch,  Phys. Rev. Lett. 122 (2019) 135001.
  * [3] SHUMLAK, U. et al., Advanced physics calculations using a multi-fluid plasma model,  Computer Phys. Comm. 182 (2011) 1767–1770.
  * [4] LEVITT, B. et al., The Zap Energy approach to commercial fusion,  Phys. Plasmas 30 (in press).
  * [5] MEIER, E. T. and SHUMLAK, U., Development of five-moment two-fluid modeling for Z-pinch physics,  Phys. Plasmas 28 (2021) 092512.
  * [6] TOTH, G. et al., Scaling the ion inertial length and its implications for modeling reconnection in global simulations,  J. Geophys. Res. Space Phys. 122 (2017) 10,336–10,355.
  * [7] DATTA, I. A. M. and SHUMLAK, U., Computationally efficient high-fidelity plasma simulations by coupling multi-species kinetic and multi-fluid models on decomposed domains,  J. Comput. Phys. 483 (2023) 112073.
  * [8] MIKELLIDES, I. G., TASSIS, K., and YORKE, H. W., 2D magnetohydrodynamic simulations of induced plasma dynamics in the near-core region of a galaxy cluster,  Mon. Not. R. Astron. Soc. 410 (2011) 2602–2616.
  * [9] VONNEUMANN, J. and RICHTMYER, R. D., A Method for the Numerical Calculation of Hydrodynamic Shocks,  J. Appl. Phys. 21 (1950) 232–237.
  * [10] SHUMLAK, U., MEIER, E. T., and LEVITT, B. J., Fusion gain an triple product for the sheared-flow-stabilized Z pinch,  Fusion Sci. Technol. 0 (2023) 1–16.
  * [11] STEPANOV, A. D. et al., Flow Z-pinch plasma production on the FuZE experiment,  Phys. Plasmas 27 (2020) 112503.
  * [12] JAHN, R. G.,  Physics of Electric Propulsion,  Dover Publications, Inc., Mineola, New York, USA (2006).
  * [13] CASSIBRY, J. T., THIO, Y. C. F., and WU, S. T., Two-dimensional axisymmetric magnetohydrodynamic analysis of blow-by in a coaxial plasma accelerator,  Phys. Plasmas 13 (2006) 053101.
  * [14] BOSCH, H.-S. and HALE, G. M., Improved formulas for fusion cross-sections and thermal reactivities,  Nucl. Fusion 32 (1992) 611–631.
  * [15] SHUMLAK, U. and HARTMAN, C. W., Sheared flow stabilization of the m = 1 kink mode in Z pinches,  Phys. Rev. Lett. 75 (1995) 3285–3288.
  * [16] BRAGINSKII, S. I.,  Review of Plasma Physics, volume 1,  Consultants Bureau, New York, USA (1965).
  * [17] ROGNLIEN, T. D. and RENSINK, M. E., Edge-plasma models and characteristics for magnetic fusion energy devices,  Fusion Eng. Design 60 (2002) 497–514.
  * [18] CLAVEAU, E. L. et al., Plasma exhaust in a sheared-flow-stabilized Z pinch,  Phys. Plasmas 27 (2020) 092510.
  * [19] STANGEBY, P. C.,  The Plasma Boundary of Magnetic Fusion Devices,  Institute of Physics Publishing, Philadelphia, Pennsylvania, USA (2000).
  * [20] SHUMLAK, U. and LOVERICH, J., Approximate Riemann solver for the two-fluid plasma model,  J. Comput. Phys. 187 (2003) 620–638.
  * [21] REDDELL, N.,  A Kinetic Vlasov Model for Plasma Simulation Using Discontinuous Galerkin Method on Many-Core Architectures,  PhD thesis, University of Washington, Seattle, Washington, USA (2016).
  * [22] HO, A., DATTA, I. A. M., and SHUMLAK, U., Physics-based-adaptive plasma model for high-fidelity numerical simulations,  Front. Phys. 6 (2018) 105.
  * [23] MEIER, E. T. and SHUMLAK, U., A general nonlinear fluid model for reacting plasma-neutral mixtures,  Phys. Plasmas 19 (2012) 072508.


# Extending the fair sampling assumption using causal diagrams

Valentin Gebhart  QSTAR, INO-CNR and LENS, Largo Enrico Fermi 2, 50125
Firenze, Italy<EMAIL_ADDRESS>Augusto Smerzi  QSTAR, INO-CNR and LENS,
Largo Enrico Fermi 2, 50125 Firenze, Italy<EMAIL_ADDRESS>
###### Abstract

Discarding undesirable measurement results in Bell experiments opens the
detection loophole that prevents a conclusive demonstration of nonlocality. As
closing the detection loophole represents a major technical challenge for many
practical Bell experiments, it is customary to assume the so-called fair
sampling assumption (FSA) that, in its original form, states that the
collectively postselected statistics are a fair sample of the ideal
statistics. Here, we analyze the FSA from the viewpoint of causal inference:
We derive a causal structure that must be present in any causal model that
faithfully encapsulates the FSA. This provides an easy, intuitive, and
unifying approach that includes different accepted forms of the FSA and
underlines what is really assumed when using the FSA. We then show that the
FSA can not only be applied in scenarios with non-ideal detectors or
transmission losses, but also in ideal experiments where only parts of the
correlations are postselected, e.g., when the particles’ destinations are in a
superposition state. Finally, we demonstrate that the FSA is also applicable
in multipartite scenarios that test for (genuine) multipartite nonlocality.

##  1 Introduction

Bell nonlocality [1, 2] represents one of the central pillars of modern
research in quantum foundations and the development of quantum technologies
[3]. A widely-used technique in Bell experiments is the discard of events of
an incomplete detection such as, e.g., the non-detection of parts of the
system due to particle losses. By the selection bias [4], the postselection
opens the detection loophole, i.e., the possibility of a local hidden variable
(LHV) model to describe the observed correlations even if the postselected
correlations violate a Bell inequality [5, 6]. The detection loophole is not
just conspiratorial: It was used in experiments to create fake violations of
Bell inequalities [7, 8, 9, 10].

Ideally, one can close the detection loophole by including the non-detection
events [6, 11, 12, 13], or by sharpening the Bell inequalities [14, 15, 38].
These methods require high detection efficiencies and have recently been
implemented in sophisticated Bell experiments that close the detection
loophole [16, 17, 18] (also while simultaneously closing the locality loophole
[19, 20, 21]). However, the required detection efficiencies represent a severe
technical challenge for the practicality of Bell experiments. Thus, a widely-
used way out is to rely on the fair sampling assumption (FSA) [22, 6, 23, 24]
that is commonly known as the assumption that the postselected statistics is a
fair representation of (i.e., is identical to) the statistics that would have
been observed using perfect detectors and no losses. An alternative form of
the FSA is the assumption that the detector settings have no influence on the
detection probability of the particles. In the latter case, the postselected
statistics need not be identical to the statistics that would be observed with
ideal detectors but, nonetheless, the postselection cannot create any fake
nonlocal correlations. We emphasize that the FSA does not close the detection
loophole, but it rather represents an assumption that restricts the possible
LHV descriptions for the measurement data. Due to the high technical
requirements of closing the detection loophole (e.g., highly efficient
detectors and minimal transmission losses), the FSA is still widely used in
Bell experiments [25, 26, 27, 28, 29, 30, 31].

In this work, we analyze the FSA from the viewpoint of causal inference and
causal diagrams [4]. In particular, we ask what structure any causal diagram
of a Bell experiment must possess to allow for a valid demonstration of
nonlocality if the data is collectively postselected. This structure should
guarantee that the postselection cannot create fake nonlocal correlations due
to the selection bias. Importantly, we ensure that the causal diagram provides
a meaningful description of the experiment by disallowing any kind of fine-
tuning of causal influences [4, 32, 33], in contrast to previous studies of
the FSA. This results in an easy and intuitive way to understand different
forms of the FSA found in the literature. Our analysis highlights what is
really assumed when using the FSA, and allows us to identify Bell scenarios
where no faithful causal explanation of the FSA exists. Furthermore, we show
that the such-obtained causal-diagram FSA can be applied to different
experiments where the correlations must be postselected even in the ideal
noiseless setup, and not just in the standard setup with non-ideal detectors
and particle losses. Finally, we prove that the FSA is also useful in
experiments that demonstrate multipartite nonlocality and genuine multipartite
nonlocality.

##  2 Fair sampling in the bipartite scenario

In the following, we derive a necessary causal structure for any faithful
causal model of a bipartite Bell scenario that includes a collective
postselection, without potentially creating fake correlations that violate the
Bell inequality. The central promise is to have a causal description that does
not employ any fine-tuning of its causal influences. By a fine-tuning of the
parameters of a causal model, two variables can be made statistically
independent even though they seem to affect each other in the causal diagram.
Instead, without fine-tuning, any statistical independence between two
variables must be evident from the diagram’s structure. If one allows for
fine-tuning, the description in terms of causal diagrams becomes irrelevant
[4, 32, 33] because any statistical (in)dependence can just be directly
imposed by hand. Causal models that are not fine-tuned are usually called
faithful models.

The central tool to infer independencies from a causal diagram are the
$d$-separation rules [4] that dictate whether a given path connecting two
variables of the causal diagram is open (i.e., the variables are generally
dependent) or blocked (the variables must be independent), also when
conditioning on other variables of the diagram. The rules state that (i) a
path is blocked if along it there is a collider (a variable where two causal
arrows collide), (ii) a path is blocked if along it there is a non-collider
that is conditioned on, and (iii) a path is open if along it there is a
collider and we condition on the collider or its descendants. The last rule
manifests itself in the selection bias and the Berkson paradox [4]. We want to
note that in the context of Bell experiments, causal diagrams and the
$d$-separation rules have been used to show that quantum violations of Bell
inequalities require fine-tuning in classical-causal explanations (e.g.,
superdeterminism, superluminal influences, or retrocausal influences) [33, 34,
35], and to show that certain collective postselection strategies are safe for
the demonstration of multipartite nonlocality [36] and genuine multipartite
nonlocality [37, 38]. In the latter works, it was shown that if the
postselection can be decided by excluding some of the parties, the detection
loophole can be closed, so one does not have to rely on the FSA. In contrast,
here we consider a postselection that must be decided by all parties together,
such that the results of Refs. [36, 37, 38] do not apply.

$\Lambda$$A$$B$$X$$Y$ (a)

$\Lambda$$A$$B$$X$$Y$$K$(b)

Figure 1: (a) The causal diagram of the local hidden variable (LHV) model in
the standard bipartite Bell scenario. (b) A collective postselection is
indicated as the variable $K$ that is influenced by both Alice’s outcome $A$
and Bob’s outcome $B$, opening the detection loophole by the selection bias
[5]. The conditioning on the decision of the postselection is indicated as a
box around $K$.

In the standard bipartite Bell scenario, two parties, Alice and Bob, share the
two parts of a quantum system and each perform local measurements on their
subsystem. Alice (Bob) chooses the measurement setting $x$ ($y$) and records
the measurement outcome $a$ ($b$). The measured correlations are called local
if they can be described by a LHV model of the form [1, 2]

$p_{ab|xy}=\sum_{\lambda}p_{\lambda}p_{a|x\lambda}p_{b|y\lambda},$ (1)

where $p_{\lambda}$, $p_{a|x\lambda}$, and $p_{b|y\lambda}$ are (conditional)
probabilities, each summing to $1$, e.g., $\sum_{\lambda}p_{\lambda}=1$. The
corresponding causal diagram is shown in Fig. 1(a). Using Eq. (1), one can
derive Bell inequalities, a violation of which proves that the correlations
are nonlocal.

In Fig. 1(b), we include the variable $K$ representing the decision of the
collective postselection (e.g., $K=1$ for postselecting the results and $K=0$
for discarding the results). If the postselected statistics $p_{ab|xy1}$ ($1$
is the value of the variable $K$) can be described by a LHV similar to Eq.
(1), they must also fullfill the Bell inequality and the postselection is
valid. By the definition of conditional probability, we can write
$p_{ab|xyk}=\sum_{\lambda}p_{\lambda|xyk}p_{ab|xy\lambda k}$, such that we can
identify the two conditions

$\displaystyle\mathrm{\textbf{CI}:}$ $\displaystyle\quad\quad\quad
p_{\lambda|xyk}=p_{\lambda|k}$ $\displaystyle\mathrm{\textbf{CII}:}$
$\displaystyle\quad\quad p_{ab|xy\lambda k}=p_{a|x\lambda k}p_{b|y\lambda k}$

that ensure a LHV description of $p_{ab|xyk}$ and thus a valid
postselection111It would actually suffice to show that $p_{ab|xyk}$ admits a
LHV model only for $k=1$, and not necessarily for all $k$. The causal-
inference techniques that we employ below yield the conditions for all $k$,
implying the case $k=1$.. We note that CI and CII correspond to the
measurement-independence and locality assumptions of Bell’s theorem,
respectively. Now consider a postselection that is collective: Both parties
must consult each other to decide the postselection. Note that if each party
can decide the postselection locally, there is no need for a FSA because the
postselection is known to be safe [13, 36]. Therefore, in the causal diagram,
$K$ is influenced by both measurement outcomes $A$ and $B$ 222Strictly
speaking, a collective postselection could also include direct influences from
the settings $X$ and $Y$ to $K$, e.g., a postselection influenced only by $X$
and $B$, or only by $X$ and $Y$. A postselection influenced by $X$ and $B$
leads to a conflict with condition CI, similar a postselection influenced by
$A$ and $B$ as in the main text. On the other side, a postselection that is
decided only by $X$ and $Y$ is a safe postselection, i.e., the conditions CI
and CII still hold. This can be seen by noting that, if $K$ is influenced only
by $X$ and $Y$, one simply has $p_{ab|xyk}=p_{ab|xy}$., see Fig. 1(b).
However, a postselection described by the causal diagram in Fig. 1(b) is in
conflict with condition CI: There is an open path $X\rightarrow A\rightarrow
K\leftarrow B\leftarrow\Lambda$ ($K$ is a collider that is conditioned on),
which contradicts $p_{\lambda|xyk}=p_{\lambda|k}$.

The causal diagram of Fig. 1(b) can thus not give a causal account of the fair
sampling without employing a fine-tuning of causal influences. In the
following, we show that any causal description of the FSA requires a certain
type of structure in the causal model if the model is not fine-tuned. We
consider a general bipartite Bell scenario where Alice’s (Bob’s) measurement
settings are given by a number of setting variables
$\mathbf{x}=(x_{1},\dots,x_{n_{A}})$ [$\mathbf{y}=(y_{1},\dots,y_{n_{B}})$],
and Alice (Bob) observes a number of measurement outcomes
$\mathbf{a}=(a_{1},\dots,a_{m_{A}})$ [$\mathbf{b}=(b_{1},\dots,b_{m_{B}})$].
Their measurement outcomes are correlated via the LHV $\Lambda$ (note that we
group all LHVs into the single LHV $\Lambda$ without loss of generality). At
this point, the different outcome variables of each party can be arbitrarily
connected by causal influences, e.g., one could have a causal influence
$A_{1}\rightarrow A_{2}$. After performing their measurements, the parties
collectively postselect their data, represented as the binary postselection
variable $K$ as above. The corresponding causal structure is identical to the
one of Fig. 1(b), expect that all setting and outcome variables are replaced
by multi-variable versions. Thus, similar to above, this general causal model
is in conflict with the conditions CI and CII for a valid postselection.

$\Lambda$$\mathbf{A}^{(K)}$$\mathbf{A}^{(X)}$$\mathbf{B}^{(K)}$$\mathbf{B}^{(Y)}$$\mathbf{X}$$\mathbf{Y}$$K$(a)

$\Lambda$$\mathbf{A}^{(K)}$$\mathbf{A}^{(X)}$$\mathbf{B}^{(K)}$$\mathbf{B}^{(Y)}$$\mathbf{X}$$\mathbf{Y}$$K$(b)

$\Lambda$$\mathbf{A}^{(K)}$$\mathbf{A}^{(X)}$$\mathbf{B}^{(K)}$$\mathbf{B}^{(Y)}$$\mathbf{X}$$\mathbf{Y}$$K$(c)

Figure 2: Causal diagrams to derive the causal structure to faithfully account
for the fair sampling assumption (FSA). (a) We divide the measurement outcome
variables $\mathbf{A}$ and $\mathbf{B}$ of a general bipartite Bell scenario
into those that are used in the postselection decision ($\mathbf{A}^{(K)}$ and
$\mathbf{B}^{(K)}$) and those that are not ($\mathbf{A}^{(X)}$ and
$\mathbf{B}^{(Y)}$). Bidirected arrows indicate arbitrary causal influences
between variables. (b) Causal diagram [restricting the one in (a)] that gives
a faithful account of condition CI. (c) Causal diagram [restricting the one in
(b)] that gives a faithful account of both conditions CI and CII,
corresponding to a valid postselection and thus to the FSA.

To derive the causal structure required to show conditions CI and CII without
fine-tuning, we first divide Alice’s outcomes $\mathbf{A}$ into the outcomes
$\mathbf{A}^{(K)}$ that are used to decide the postselection $K$, and the
outcomes $\mathbf{A}^{(X)}$ that are not, see Fig. 2(a). Similarly, we divide
$\mathbf{B}$ into $\mathbf{B}^{(K)}$ and $\mathbf{B}^{(Y)}$. The bidirected
arrows depict arbitrary causal influences, i.e., $A_{i}\multimapboth A_{j}$
includes $A_{i}\rightarrow A_{j}$, $A_{i}\leftarrow A_{j}$, or a hidden
variable $\Gamma$ such that $A_{i}\leftarrow\Gamma\rightarrow A_{j}$ (a hidden
common cause), and combinations of thereof. Note that a hidden common cause
can be included in the LHV $\Lambda$. Now, if some $A_{j}\in\mathbf{A}^{(K)}$
is directly influenced by $\mathbf{X}$, there is an open path
$\mathbf{X}\rightarrow A_{j}\leftarrow\Lambda$ because we condition on $K$, a
descendant of the collider $A_{j}$ (we assume that the LHV $\Lambda$
influences all measurement outcomes). Thus, to preserve condition CI,
$\mathbf{X}$ cannot have a direct influence on the group $\mathbf{A}^{(K)}$,
and similarly for $\mathbf{Y}$ and $\mathbf{B}^{(K)}$. Next, if there was an
influence $A_{i}\rightarrow A_{j}$ with $A_{i}\in\mathbf{A}^{(X)}$ and
$A_{j}\in\mathbf{A}^{(K)}$, there would be an open path $X\rightarrow
A_{i}\rightarrow A_{j}\leftarrow\Lambda$, in conflict with CI. Here we assume
that $\mathbf{X}$ influences $A_{i}$ because, otherwise, $A_{i}$ would neither
be useful for violating Bell inequalities333 To violate Bell inequalities,
each party’s setting must influence its measurement outcome. For instance, in
the bipartite scenario with one setting and one measurement variable per
party, assume that Alice’s setting does not influence her outcome. Due to non-
signalling, her setting cannot influence Bob’s outcome, so one has
$p_{ab|xy}=p_{ab|y}=p_{a|y}p_{b|ay}=p_{a}p_{b|ay},$ where we have used the no-
signalling principle $p_{a|y}=p_{a}$. This yields a LHV model for $p_{ab|xy}$:
Defining $\Lambda$ to take the same values as $A$ with identical
probabilities, $p_{\lambda}=p_{a}$ for $\lambda=a$, one has
$p_{ab|xy}=p_{a}p_{b|ay}=\sum_{\lambda}p_{\lambda}\delta_{a,\lambda}p_{b|\lambda,y},$
where $\delta$ is the Kronecker symbol. , nor would it be useful to decide the
postselection (because of $A_{i}\nrightarrow K$), so one can just discard the
outcome $A_{i}$ from the analysis. We thus obtain Fig. 2(b) which ensures that
the condition CI is fulfilled. For instance, for proving that
$p_{\lambda|\mathbf{x}k}=p_{\lambda|k}$, note that the only path connecting
$\mathbf{X}$ and $\Lambda$ passes through $\mathbf{A}^{(X)}$ that, being a
collider, blocks the path.

Since there is no influence from $X$ to any $A_{j}\in\mathbf{A}^{(K)}$, the
variables in $\mathbf{A}^{(K)}$ (and $\mathbf{B}^{(K)}$) are not useful to
violate a Bell inequality${}^{\textrm{\ref{footnote3}}}$. We thus only
consider the variables $\mathbf{A}^{(X)}$ and $\mathbf{B}^{(Y)}$ as inputs to
the Bell inequality. To show that the postselected statistics (describing
correlations between $\mathbf{A}^{(X)}$ and $\mathbf{B}^{(Y)}$) can be
described by a LHV model, it remains to show condition CII, i.e.,
$p_{a_{i}b_{r}|\mathbf{x}\mathbf{y}\lambda k}=p_{a_{i}|\mathbf{x}\lambda
k}p_{b_{r}|\mathbf{y}\lambda k}$ for all $A_{i}\in\mathbf{A}^{(X)}$ and
$B_{r}\in\mathbf{B}^{(Y)}$. If there were influences $A_{j}\rightarrow A_{i}$
and $B_{s}\rightarrow B_{r}$ for $A_{i}\in\mathbf{A}^{(X)}$,
$A_{j}\in\mathbf{A}^{(K)}$, $B_{r}\in\mathbf{B}^{(Y)}$ and
$B_{s}\in\mathbf{B}^{(K)}$, the path $A_{i}\leftarrow A_{j}\rightarrow
K\leftarrow B_{s}\rightarrow B_{r}$ would be open, in conflict with CII.
Excluding influences of the form $\mathbf{A}^{(K)}\rightarrow\mathbf{A}^{(X)}$
and $\mathbf{B}^{(K)}\rightarrow\mathbf{B}^{(Y)}$ ensures CII because the only
paths that connect $\mathbf{A}^{(X)}$ and $\mathbf{X}$ to $\mathbf{B}^{(Y)}$
and $\mathbf{Y}$ pass though $\Lambda$ and are blocked because $\Lambda$ is a
non-collider that is conditioned on. Thus, we conclude with Fig. 2(c) which
ensures both CI and CII (if restricted to variables in $\mathbf{A}^{(X)}$ and
$\mathbf{B}^{(Y)}$).

In summary, we have started from a general bipartite Bell scenario including a
collective postselection and derived a necessary causal structure to
faithfully describe the FSA, see Fig. 2(c). This structure requires that each
party must have at least one measurement variable ($\mathbf{A}^{(K)}$ and
$\mathbf{B}^{(K)}$) that is used to decide the postselection and that is
independent of the measurement settings, and at least one measurement variable
($\mathbf{A}^{(X)}$ and $\mathbf{B}^{(Y)}$) that is used as an input in the
Bell inequality and that does not influence the postselection. The smallest
realization of this structure is a Bell scenario where each party has a binary
measurement setting ($X$ and $Y$), a binary measurement variable ($A_{2}$ and
$B_{2}$) that dictates the postselection and a binary measurement variable
($A_{1}$ and $B_{1}$) that is used in the Bell inequality. In the standard use
of the FSA to deal with the loss of particles, $A_{2}$ and $B_{2}$ correspond
to the number of detected particles in the respective measurement stations,
while the $A_{1}$ and $B_{1}$ correspond to the outcomes of, e.g., a
polarization measurement of incoming photons444 In this example, one could
wonder about the meaning of $A_{1}$ if Alice does not detect a particle
($A_{2}=0$). Since this event will be discarded in the postselection, the
value attributed to $A_{1}$ is not important. To be consistent with the
assumption of Fig. 2(c) ($A_{2}\nrightarrow A_{1}$) one could, e.g., flip a
coin to set the value of $A_{1}$ in this case.. We note however that, in
general, the variables $\mathbf{A}^{(K)}$ and $\mathbf{B}^{(K)}$ not
necessarily represent the number of detected particles, but the results hold
for any collectively decided postselection.

We want to emphasize that the causal-diagram FSA is not only applicable to the
standard scenario where one particle is sent to each party but there are
detection and transmission losses, but also for certain postselection methods
if the particles are generated in a superposition of their destinations [13,
39, 40]. In particular, demonstrations that a coincidence postselection is
safe if the number of particles is conserved [36, 37, 38] become unnecessary
if one assumes the FSA. In other words, the FSA covers both a postselection
due to inefficient detectors and transmission losses, and a postselection in
ideal experiments due to a varying distribution of particles.

Above, we have derived a necessary structure of any bipartite causal diagram
that faithfully describes the FSA for a general collective postselection.
However, typical applications of the FSA are situations in which each
measurement party needs to detect a single particle. Here, the variables
$A_{2}$ and $B_{2}$ correspond to the number of detected particles in Alice’s
and Bob’s experiment, respectively. This is a special case of a collective
postselection, in which postselecting an event (denoted as $K=1$ above) is
equivalent to a fixed combination of values for $A_{2}$ and $B_{2}$
($A_{2}=B_{2}=1$). Thus, one can simply use a causal diagram with a
conditioning on the variables $A_{2}$ and $B_{2}$ without introducing the
postselection variable $K$. While influences such as $X\rightarrow A_{2}$ and
$A_{1}\rightarrow A_{2}$ are still in conflict with condition CI, an influence
of the form $A_{2}\rightarrow A_{1}$ can now be allowed. The corresponding
causal diagram is shown in Fig. 3. The causal diagram of Fig. 2(c) (and its
smallest realization) is more general though: An example of a collective
postselection that cannot be modeled with Fig. 3 is when the parties
postselect events for which $A_{2}=B_{2}$. Here, a postselected event does not
imply fixed values of $A_{2}$ and $B_{2}$.

$\Lambda$$A_{2}$$A_{1}$$B_{2}$$B_{1}$$X$$Y$ Figure 3: A possible causal
diagram for the FSA if, for the postselection of the results, each party must
receive a single particle. Here, $A_{2}$ ($B_{2}$) corresponds to the number
of particles detected in Alice’s (Bob’s) measurement device.

To conclude this section, we want to mention that, while one cannot
experimentally certify that the LHV is causally modelled by Fig. 2(c) (in the
same way that one cannot experimentally certify the FSA), our results
highlight in which cases is it not possible to have a faithful causal account
of the FSA. As discussed above, to account for the FSA, each party must have
at least two separate measurement results, one that influences the
postselection decision, and one that is used in the Bell inequality. This
excludes setups for which the outcome used in the Bell inequality is also used
to decide the postselection. An example is a setup proposed by Franson [41] to
create nonlocality from energy-time entanglement, which has been shown to
admit LHV models that reproduce the observed statistics and the apparent Bell
inequality violation even in the noiseless case [42, 43, 44, 45]. Here, each
party has two different measurement outcomes, an early detection time and a
late detection time. Even without particle losses, the statistics must be
collectively postselected in order to violate a Bell inequality: Only those
events are postselected for which both particles arrive either at the early or
at the late detection time. Since the time of arrival is used both in the
postselection and in the Bell inequality, there is no way to introduce two
separate variables per party, e.g., $A_{1}$ and $A_{2}$, with the roles as
above. Thus, a FSA for the noiseless version of the original Franson setup
must rely on a fine-tuning in the causal diagram.

###  2.1 Comparison to standard FSAs

We now want to briefly compare the causal-diagram FSA to its different forms
found in the literature. We emphasize that any of the following forms of the
FSA corresponds to a fine-tuning condition on the original causal diagram of
the Bell experiment (Fig. 1). We first comment on the common
(mis-)understanding that the FSA means that there is no observable influence
of the measurement setting on the probability of detecting a particle,
$p_{d|x}=p_{d}$, where $d$ represents Alice’s detection of a particle. As
shown in Ref. [23] with a counter example, this assumption does not ensure a
safe postselection. The original way of stating the FSA is that the
postselected statistics should be a fair sample of (i.e., be identical to) the
statistics that would have been obtained using perfect detectors [6]. This
assumption is satisfied if $p_{d|x\lambda}=p_{d}$, i.e., if the probability of
detecting a particle depends neither on the setting $X$ nor on the LHV
$\Lambda$. This condition ensures a safe postselection but, as we have seen
above, can be weakened: Assuming that $p_{d|x\lambda}=p_{d|\lambda}$ already
provides a safe postselection [23]. Here, the postselected ensemble may differ
from the original one, $p_{\lambda|k}\neq p_{\lambda}$, but it still has the
form of a LHV model, Eq. (1). Assuming a causal-diagram representation without
fine-tuning, the assumption that $X$ cannot influence the detection variable
corresponds to the above causal diagrams of Fig. 2(c) or Fig. 3.

Finally, we note that this FSA can be further weakened to the assumption that
$p_{d|x\lambda}=\eta^{(d)}_{x}\eta^{(d)}_{\lambda}$ [23, 24], i.e., the
assumption that the detection efficiency depends on both the setting $X$ and
the LHV $\Lambda$ but it factorizes. This factorization condition cannot be
depicted in a causal diagram and, as the dependence implies that $A_{2}$ is
influenced by both $X$ and $\Lambda$, it represents a fine-tuning of the
causal influences.

##  3 Fair sampling for genuine multipartite nonlocality

The causal-diagram FSA of Fig. 2(c) (or Fig. 3) is also applicable in
multipartite Bell experiments. For more than two measurement parties, there
are different notions of nonlocality that can be demonstrated by a violation
of the corresponding inequalities [46, 47, 48, 49, 50]. For simplicity, we
focus on the three-partite case, but the discussion holds for any number of
parties. We thus include a third party Charlie who chooses a measurement
setting $z$ and observes the measurement outcomes $(c_{1},c_{2})$, where for
simplicity we only consider the smallest realization of the FSA causal
structure including two measurement variables per party. First, one can assume
a LHV model in the multipartite case similar to Eq. (1), corresponding to a
causal structure as shown in Fig. 4(a), where we included the FSA derived
above. Using the LHV model, one can demonstrate inequalities that test
multipartite nonlocality [47]. The validity of the postselection, namely the
conditions $p_{\lambda|xyzk}=p_{\lambda|k}$ and $p_{a_{1}b_{1}c_{1}|xyz\lambda
k}=p_{a_{1}|x\lambda k}p_{b_{1}|y\lambda k}p_{c_{1}|z\lambda k}$, can be shown
in exact analogy to the bipartite case above.

$X$(a)$A_{1}$$A_{2}$$\Lambda$$B_{2}$$B_{1}$$Y$$C_{1}$$C_{2}$$Z$$K$

$X$(b)$A_{1}$$A_{2}$$\Lambda$$B_{2}$$B_{1}$$Y$$C_{1}$$C_{2}$$Z$$K$

$X$(c)$A_{1}$$A_{2}$$\Lambda$$B_{2}$$B_{1}$$Y$$C_{1}$$C_{2}$$Z$$K$

Figure 4: Examples of the FSA causal structure in the three-partite Bell
scenario, for (a) a local hidden variable (LHV) model that corresponds to
tests of multipartite nonlocality and (b,c) a hybrid local-nonlocal hidden
variable model that corresponds to tests of genuine multipartite nonlocality.
(c) In the hybrid model, when conditioning on a specific value of the LHV
$\Lambda$, two of the three parties can share nonlocal quantum correlations
(light blue lines) that are subject to the no-signalling fine-tuning
condition.

A second and stronger form of three-partite nonlocality is genuine three-
partite nonlocality. Here, instead of assuming a LHV model, one allows for two
of the three parties to share nonlocal quantum correlations, in a model that
is called a hybrid local-nonlocal hidden variable model [46]555We note that,
recently, several new definitions of (genuine) multipartite nonlocality have
been proposed that we do not address in this work. These include network
nonlocality [51], broadcasting correlations [52], and nonlocality that is
based on the resource theory of local operations and shared randomness (LOSR)
[53, 54]. See also Ref. [55] for discussion of different classes of
multipartite nonlocality.. These quantum correlations, fulfilling the no-
signalling principle [49], cannot be depicted in a classical causal diagram
without using fine-tuning conditions [33, 34]. In Fig. 4, we indicate these
correlations as light blue lines between the outcome variables, reminding that
these influences are subject to the no-signalling principle. The hybrid model
then dictates that, given a specific hidden variable $\lambda$ (i.e., when
conditioning on $\Lambda$), there can only be nonlocal correlations between
two of the parties, see Fig. 4(c). In contrast, when not conditioning on
$\Lambda$, there can possibly exist nonlocal correlations between any pair of
parties, see Fig. 4(b), where we use different colors to emphasize that only
one pair of the parties can share nonlocal correlations at a time.

Hybrid local-nonlocal hidden variable models fulfill certain inequalities that
test for genuine multipartite nonlocality [46], and, similar to above, there
are conditions on the postselected statistics that, if fulfilled, prove that a
collective postselection is safe [37]. One can directly show that, using the
causal-diagram FSA as shown in Fig. 4, the conditions of a safe postselection
are fulfilled. For instance, for the first condition,
$p_{\lambda|xyzk}=p_{\lambda|k}$, we note that a path such as $X\rightarrow
A_{1}\rightarrow B_{2}\rightarrow K\leftarrow C_{2}\leftarrow\Lambda$, that
appears to be an open path since the collider $K$ is conditioned on, is
blocked due to the no-signalling condition: Alice’s measurement setting $X$
cannot influence Bob’s measurement outcome $B_{2}$. The validity of the second
condition of a factorization given a specific value of $\Lambda$, e.g.,
$p_{a_{1}b_{1}c_{1}|xyz\lambda k}=p_{a_{1}b_{1}|xy\lambda k}p_{c_{1}|z\lambda
k}$, see Fig. 4(c), can again be seen by noting that any path that connects
$C_{1}$ and $Z$ to the other parties passes through $\Lambda$, and is thus
blocked because $\Lambda$ is a non-collider that is conditioned on. Thus we
see that the FSA depicted in Fig. 4 also suffices to validate a collective
postselection for the demonstration of genuine multipartite nonlocality.

Finally, as in the bipartite case, if each party needs to measure a single
particle, the FSA can also be explained by a three-partite causal diagram
similar to Fig. 3. Here, $A_{2}$ represents the number detected particles in
Alice’s measurement, and one can allow for an influence of the form
$A_{2}\rightarrow A_{1}$, and similarly for Bob and Charlie. Also this version
of a fair-sampling causal diagram suffices to prove the conditions of a safe
postselection for both multipartite nonlocality and genuine multipartite
nonlocality.

##  4 Conclusions

We have discussed a causal explanation of the fair sampling assumption (FSA)
that ensures that a collective postselection cannot create fake nonlocal
correlations in Bell experiments. For this purpose, we have derived a causal
structure that any causal model of a bipartite Bell scenario must possess to
guarantee that the postselected statistics take the form of a local hidden
variable (LHV) model, without requiring a fine-tuning of the causal
influences. We have employed the framework of causal inference and
$d$-separation rules [4] as a mediator between a causal structure and the
implied relations of conditional independence between the variables. Our
results clarify what one really assumes if one uses the FSA, and while the
corresponding causal structure is not experimentally certifiable (similar to
the FSA itself), our results demonstrate that, in certain Bell scenarios,
there is no faithful causal account of the FSA. The derived FSA causal
structure yields an easy and intuitive explanation of the FSA and can be used
to understand different forms of the FSA in the literature. Furthermore,
besides standard Bell scenarios with a non-ideal detection and particle
losses, we have demonstrated that the causal-diagram FSA can also be applied
in noise-free scenarios where the statistics must be postselected because the
particles are randomly distributed between the parties [39, 40, 37, 38].
Finally, we have shown that the FSA is also applicable for a collective
postselection in demonstrations of multipartite nonlocality and genuine
multipartite nonlocality.

## Acknowledgments

This work was supported by the European Commission through the H2020 QuantERA
ERA-NET Cofund in Quantum Technologies project “MENTA”

## References

  * [1] John S Bell.  “On the einstein podolsky rosen paradox”.  Physics 1, 195 (1964).
  * [2] John S Bell.  “The theory of local beables”.  In Speakable and Unspeakable in Quantum Mechanics: Collected Papers on Quantum Philosophy.  Pages 52–62.  Cambridge University Press (2004).  2 edition.
  * [3] Nicolas Brunner, Daniel Cavalcanti, Stefano Pironio, Valerio Scarani, and Stephanie Wehner.  “Bell nonlocality”.  Rev. Mod. Phys. 86, 419–478 (2014).
  * [4] Judea Pearl.  “Causality: Models, reasoning, and inference”.  Cambridge University Press. (2009).
  * [5] Philip M. Pearle.  “Hidden-variable example based upon data rejection”.  Phys. Rev. D 2, 1418–1425 (1970).
  * [6] John F. Clauser and Michael A. Horne.  “Experimental consequences of objective local theories”.  Phys. Rev. D 10, 526–535 (1974).
  * [7] D. S. Tasca, S. P. Walborn, F. Toscano, and P. H. Souto Ribeiro.  “Observation of tunable popescu-rohrlich correlations through postselection of a gaussian state”.  Phys. Rev. A 80, 030101 (2009).
  * [8] Ilja Gerhardt, Qin Liu, Antía Lamas-Linares, Johannes Skaar, Valerio Scarani, Vadim Makarov, and Christian Kurtsiefer.  “Experimentally faking the violation of bell’s inequalities”.  Phys. Rev. Lett. 107, 170404 (2011).
  * [9] Enrico Pomarico, Bruno Sanguinetti, Pavel Sekatski, Hugo Zbinden, and Nicolas Gisin.  “Experimental amplification of an entangled photon: what if the detection loophole is ignored?”.  New J. Phys. 13, 063031 (2011).
  * [10] J Romero, D Giovannini, D S Tasca, S M Barnett, and M J Padgett.  “Tailored two-photon correlation and fair-sampling: a cautionary tale”.  New J. Phys. 15, 083047 (2013).
  * [11] N. David Mermin.  “The EPR experiment—thoughts about the “loophole””.  Ann. N. Y. Acad. Sci. 480, 422–427 (1986).
  * [12] Philippe H. Eberhard.  “Background level and counter efficiencies required for a loophole-free einstein-podolsky-rosen experiment”.  Phys. Rev. A 47, R747–R750 (1993).
  * [13] Fabio Sciarrino, Giuseppe Vallone, Adán Cabello, and Paolo Mataloni.  “Bell experiments with random destination sources”.  Phys. Rev. A 83, 032112 (2011).
  * [14] Anupam Garg and N. D. Mermin.  “Detector inefficiencies in the einstein-podolsky-rosen experiment”.  Phys. Rev. D 35, 3831–3835 (1987).
  * [15] Jan-Åke Larsson.  “Bell’s inequality and detector inefficiency”.  Phys. Rev. A 57, 3304–3308 (1998).
  * [16] Mary A Rowe, David Kielpinski, Volker Meyer, Charles A Sackett, Wayne M Itano, Christopher Monroe, and David J Wineland.  “Experimental violation of a bell’s inequality with efficient detection”.  Nature 409, 791–794 (2001).
  * [17] D. N. Matsukevich, P. Maunz, D. L. Moehring, S. Olmschenk, and C. Monroe.  “Bell inequality violation with two remote atomic qubits”.  Phys. Rev. Lett. 100, 150404 (2008).
  * [18] B. G. Christensen, K. T. McCusker, J. B. Altepeter, B. Calkins, T. Gerrits, A. E. Lita, A. Miller, L. K. Shalm, Y. Zhang, S. W. Nam, N. Brunner, C. C. W. Lim, N. Gisin, and P. G. Kwiat.  “Detection-loophole-free test of quantum nonlocality, and applications”.  Phys. Rev. Lett. 111, 130406 (2013).
  * [19] Lynden K. Shalm, Evan Meyer-Scott, Bradley G. Christensen, Peter Bierhorst, Michael A. Wayne, Martin J. Stevens, Thomas Gerrits, Scott Glancy, Deny R. Hamel, Michael S. Allman, Kevin J. Coakley, Shellee D. Dyer, Carson Hodge, Adriana E. Lita, Varun B. Verma, Camilla Lambrocco, Edward Tortorici, Alan L. Migdall, Yanbao Zhang, Daniel R. Kumor, William H. Farr, Francesco Marsili, Matthew D. Shaw, Jeffrey A. Stern, Carlos Abellán, Waldimar Amaya, Valerio Pruneri, Thomas Jennewein, Morgan W. Mitchell, Paul G. Kwiat, Joshua C. Bienfang, Richard P. Mirin, Emanuel Knill, and Sae Woo Nam.  “Strong loophole-free test of local realism”.  Phys. Rev. Lett. 115, 250402 (2015).
  * [20] Marissa Giustina, Marijn A. M. Versteegh, Sören Wengerowsky, Johannes Handsteiner, Armin Hochrainer, Kevin Phelan, Fabian Steinlechner, Johannes Kofler, Jan-Åke Larsson, Carlos Abellán, Waldimar Amaya, Valerio Pruneri, Morgan W. Mitchell, Jörn Beyer, Thomas Gerrits, Adriana E. Lita, Lynden K. Shalm, Sae Woo Nam, Thomas Scheidl, Rupert Ursin, Bernhard Wittmann, and Anton Zeilinger.  “Significant-loophole-free test of bell’s theorem with entangled photons”.  Phys. Rev. Lett. 115, 250401 (2015).
  * [21] Bas Hensen, Hannes Bernien, Anaïs E Dréau, Andreas Reiserer, Norbert Kalb, Machiel S Blok, Just Ruitenberg, Raymond FL Vermeulen, Raymond N Schouten, Carlos Abellán, et al.  “Loophole-free bell inequality violation using electron spins separated by 1.3 kilometres”.  Nature 526, 682–686 (2015).
  * [22] John F. Clauser, Michael A. Horne, Abner Shimony, and Richard A. Holt.  “Proposed experiment to test local hidden-variable theories”.  Phys. Rev. Lett. 23, 880–884 (1969).
  * [23] Dominic W. Berry, Hyunseok Jeong, Magdalena Stobińska, and Timothy C. Ralph.  “Fair-sampling assumption is not necessary for testing local realism”.  Phys. Rev. A 81, 012109 (2010).
  * [24] Davide Orsucci, Jean-Daniel Bancal, Nicolas Sangouard, and Pavel Sekatski.  “How post-selection affects device-independent claims under the fair sampling assumption”.  Quantum 4, 238 (2020).
  * [25] Igor Marinković, Andreas Wallucks, Ralf Riedinger, Sungkun Hong, Markus Aspelmeyer, and Simon Gröblacher.  “Optomechanical bell test”.  Phys. Rev. Lett. 121, 220404 (2018).
  * [26] Dominik Rauch, Johannes Handsteiner, Armin Hochrainer, Jason Gallicchio, Andrew S. Friedman, Calvin Leung, Bo Liu, Lukas Bulla, Sebastian Ecker, Fabian Steinlechner, Rupert Ursin, Beili Hu, David Leon, Chris Benn, Adriano Ghedina, Massimo Cecconi, Alan H. Guth, David I. Kaiser, Thomas Scheidl, and Anton Zeilinger.  “Cosmic bell test using random measurement settings from high-redshift quasars”.  Phys. Rev. Lett. 121, 080403 (2018).
  * [27] Emanuele Polino, Iris Agresti, Davide Poderini, Gonzalo Carvacho, Giorgio Milani, Gabriela Barreto Lemos, Rafael Chaves, and Fabio Sciarrino.  “Device-independent test of a delayed choice experiment”.  Phys. Rev. A 100, 022111 (2019).
  * [28] S. Gómez, A. Mattar, I. Machuca, E. S. Gómez, D. Cavalcanti, O. Jiménez Farías, A. Acín, and G. Lima.  “Experimental investigation of partially entangled states for device-independent randomness generation and self-testing protocols”.  Phys. Rev. A 99, 032108 (2019).
  * [29] Davide Poderini, Iris Agresti, Guglielmo Marchese, Emanuele Polino, Taira Giordani, Alessia Suprano, Mauro Valeri, Giorgio Milani, Nicolò Spagnolo, Gonzalo Carvacho, et al.  “Experimental violation of n-locality in a star quantum network”.  Nat. Commun. 11, 1–8 (2020).
  * [30] Santiago Tarrago Velez, Vivishek Sudhir, Nicolas Sangouard, and Christophe Galland.  “Bell correlations between light and vibration at ambient conditions”.  Sci. Adv. 6, eabb0260 (2020).
  * [31] Iris Agresti, Davide Poderini, Leonardo Guerini, Michele Mancusi, Gonzalo Carvacho, Leandro Aolita, Daniel Cavalcanti, Rafael Chaves, and Fabio Sciarrino.  “Experimental device-independent certified randomness generation with an instrumental causal structure”.  Commun. Phys. 3, 1–7 (2020).
  * [32] Peter Spirtes, Clark N Glymour, Richard Scheines, and David Heckerman.  “Causation, prediction, and search”.  MIT press. (2000).
  * [33] Christopher J Wood and Robert W Spekkens.  “The lesson of causal discovery algorithms for quantum correlations: causal explanations of bell-inequality violations require fine-tuning”.  New J. Phys. 17, 033002 (2015).
  * [34] John-Mark A. Allen, Jonathan Barrett, Dominic C. Horsman, Ciarán M. Lee, and Robert W. Spekkens.  “Quantum common causes and quantum causal models”.  Phys. Rev. X 7, 031021 (2017).
  * [35] Eric G. Cavalcanti.  “Classical causal models for bell and kochen-specker inequality violations require fine-tuning”.  Phys. Rev. X 8, 021018 (2018).
  * [36] Pawel Blasiak, Ewa Borsuk, and Marcin Markiewicz.  “On safe post-selection for Bell tests with ideal detectors: Causal diagram approach”.  Quantum 5, 575 (2021).
  * [37] Valentin Gebhart, Luca Pezzè, and Augusto Smerzi.  “Genuine multipartite nonlocality with causal-diagram postselection”.  Phys. Rev. Lett. 127, 140401 (2021).
  * [38] Valentin Gebhart and Augusto Smerzi.  “Coincidence postselection for genuine multipartite nonlocality: Causal diagrams and threshold efficiencies” (2022).  Phys. Rev. A 106, 062202 (2022).
  * [39] Bernard Yurke and David Stoler.  “Bell’s-inequality experiments using independent-particle sources”.  Phys. Rev. A 46, 2229–2234 (1992).
  * [40] Bernard Yurke and David Stoler.  “Einstein-podolsky-rosen effects from independent particle sources”.  Phys. Rev. Lett. 68, 1251–1254 (1992).
  * [41] J. D. Franson.  “Bell inequality for position and time”.  Phys. Rev. Lett. 62, 2205–2208 (1989).
  * [42] Sven Aerts, Paul Kwiat, Jan-Åke Larsson, and Marek Z˙ukowski.  “Two-photon franson-type experiments and local realism”.  Phys. Rev. Lett. 83, 2872–2875 (1999).
  * [43] Jonathan Jogenfors, Ashraf Mohamed Elhassan, Johan Ahrens, Mohamed Bourennane, and Jan Åke Larsson.  “Hacking the bell test using classical light in energy-time entanglement–based quantum key distribution”.  Sci. Adv. 1, e1500793 (2015).
  * [44] Adán Cabello, Alessandro Rossi, Giuseppe Vallone, Francesco De Martini, and Paolo Mataloni.  “Proposed bell experiment with genuine energy-time entanglement”.  Phys. Rev. Lett. 102, 040401 (2009).
  * [45] G. Lima, G. Vallone, A. Chiuri, A. Cabello, and P. Mataloni.  “Experimental bell-inequality violation without the postselection loophole”.  Phys. Rev. A 81, 040101 (2010).
  * [46] George Svetlichny.  “Distinguishing three-body from two-body nonseparability by a bell-type inequality”.  Phys. Rev. D 35, 3066–3069 (1987).
  * [47] N. David Mermin.  “Extreme quantum entanglement in a superposition of macroscopically distinct states”.  Phys. Rev. Lett. 65, 1838–1840 (1990).
  * [48] Jean-Daniel Bancal, Cyril Branciard, Nicolas Gisin, and Stefano Pironio.  “Quantifying multipartite nonlocality”.  Phys. Rev. Lett. 103, 090503 (2009).
  * [49] Jean-Daniel Bancal, Jonathan Barrett, Nicolas Gisin, and Stefano Pironio.  “Definitions of multipartite nonlocality”.  Phys. Rev. A 88, 014102 (2013).
  * [50] Patricia Contreras-Tejada, Carlos Palazuelos, and Julio I. de Vicente.  “Genuine multipartite nonlocality is intrinsic to quantum networks”.  Phys. Rev. Lett. 126, 040501 (2021).
  * [51] Miguel Navascués, Elie Wolfe, Denis Rosset, and Alejandro Pozas-Kerstjens.  “Genuine network multipartite entanglement”.  Phys. Rev. Lett. 125, 240505 (2020).
  * [52] Debashis Saha and Marcin Pawłowski.  “Structure of quantum and broadcasting nonlocal correlations”.  Phys. Rev. A 92, 062129 (2015).
  * [53] David Schmid, Thomas C. Fraser, Ravi Kunjwal, Ana Belen Sainz, Elie Wolfe, and Robert W. Spekkens.  “Understanding the interplay of entanglement and nonlocality: motivating and developing a new branch of entanglement theory” (2021).  arXiv:2004.09194.
  * [54] Xavier Coiteux-Roy, Elie Wolfe, and Marc-Olivier Renou.  “No bipartite-nonlocal causal theory can explain nature’s correlations”.  Phys. Rev. Lett. 127, 200401 (2021).
  * [55] Rafael Chaves, Daniel Cavalcanti, and Leandro Aolita.  “Causal hierarchy of multipartite Bell nonlocality”.  Quantum 1, 23 (2017).


# Pyramid Fusion Dark Channel Prior for Single Image Dehazing

###### Abstract

In this paper, we propose the pyramid fusion dark channel prior (PF-DCP) for
single image dehazing. Based on the well-known Dark Channel Prior (DCP), we
introduce an easy yet effective approach PF-DCP by employing the DCP algorithm
at a pyramid of multi-scale images to alleviate the problem of patch size
selection. In this case, we obtain the final transmission map by fusing
transmission maps at each level to recover a high-quality haze-free image.
Experiments on RESIDE SOTS show that PF-DCP not only outperforms the
traditional prior-based methods with a large margin but also achieves
comparable or even better results of state-of-art deep learning approaches.
Furthermore, the visual quality is also greatly improved with much fewer color
distortions and halo artifacts.

Index Terms—  Traditional method, single image dehazing, image pyramid, dark
channel prior

##  1 Introduction

Haze is ubiquitous in life, which comes from fog, sandstorms, or pollutants in
the air. Haze in the air largely degrades the quality of the photo. Many
details of the scene become blurry and the visibility is reduced. From the
aspect of aesthetic, the images without haze are much clearer and visually
pleasing than the hazy ones. For these reasons, haze removal has been studied
for decades, yet still a challenging task in computer vision.

The key problem of haze removal is the strategy to dehaze an image. Hautiére
et al. [1] proposed an equation to model hazy images captured by the camera:

$I(x)=J(x)t(x)+A(1-t(x)),$ (1)

where $I$ is the observed hazy image, $J$ is the real scene radiance, $A$ is
the atmospheric light, $t$ is the transmission rate. The goal of the image
dehazing is to recover haze-free scene $J$, with only the knowledge of $I$. If
we are able to estimate $A$ and $t$, we can recover the real scene $J$ as
follows:

$J(x)=\frac{I(x)-A}{t(x)}+A.$ (2)

However, $A$ and $t$ are both unknown, and it is impossible to estimate their
values from a single hazy image without any constraint. Early approaches
developed a strong prior or assumption of the difference between the haze-free
and hazy images. Tan [2] maximized the local contrast to restore the clear
images given the fact that clear images have higher contrast compared to hazy
ones. He et al. [3] utilized the prior knowledge that clear images have a dark
channel whereas hazy images do not. Meng et al. [4] introduced a boundary
constraint to help formulate the transmission map. [5] proposed a color
attenuation prior (CAP), which models a linear relationship between the scene
depth, brightness, and saturation. Berman et al. [6] suggested a non-local
prior (NLP), using only 500 distinct colors to quantize the image and
clustering different pixels according to their colors into haze-lines for
dehazing.

In the emerging years of deep learning, researchers focus on building models
to generate the transmission map, $t(x)$, and atmospheric light, $A$,
separately. Cai et al. [7] used a convolutional neural network (CNN) to
estimate the transmission rate and use the atmospheric scattering model to
dehaze the image. Ren et al. [8] introduced a multi-scale CNN model to improve
the quality of the transmission map. However, estimating two parameters
individually might cause larger errors. Li et al. [9] unified the transmission
map and the atmospheric light into one parameter and use it to estimate the
clear image. Later, Ren et al. [10] designed a coarse-to-fine fusion network,
which recovers the clear image directly from the hazy one. Deep learning has
gained huge success by surpassing the physical models built by humans and
exploring the tremendous power of computation. Nevertheless, it is troublesome
to train a deep learning model so that it performs well. Not only a large
labeled dataset is required, but also the training process is time-consuming,
let alone the problems of interpretability and robustness [11], which is
required to be employed for productization.

Among all traditional algorithms, due to the effectiveness and robustness, DCP
serves as a building block to many haze-removal tasks [12, 13]. Besides, DCP
is also broadly used as a standard baseline of the performance evaluation in
recent works using deep learning [7, 9, 10, 14, 15, 16]. However, DCP still
suffers from the problem of patch size selection. Inappropriate patch size
would introduce an inaccurate estimation of the global atmospheric light and a
poor transmission map [17], which leads to dis-satisfactory haze removal
results.

In this paper, we propose the pyramid fusion dark channel prior (PF-DCP) using
an iterative weighted fusion strategy to address the problem of patch size
selection. Through incorporating multi-scale information of the transmission
map, we demonstrate that the proposed method outperforms DCP by a large margin
with a similar computational cost. The halo artifacts of the original DCP are
greatly reduced, and the dehazed images become much smoother in colors and
visually pleasing. Although the idea of using multi-scale information is well-
established in various computer vision tasks, its application as a traditional
algorithm in this field has not been proposed yet. The novelty of this paper
lies in two aspects by solving two critical problems. First, there is no
evaluation metric to select the best patch size. Thus, we adopt pyramid fusion
among different sizes to alleviate the problem. Second, the process of
transmission map refinement helps maintain a similar level of sharpness at the
final dehazed images, which prevents the algorithm from degrading the image
quality. Furthermore, as a traditional algorithm that is interpretable and
robust, PF-DCP manages to produce superior results compared to other prior-
based traditional algorithms and commensurable to the latest deep learning
methods.

##  2 Method

###  2.1 Preliminaries

He et al. [3] proposed the Dark Channel Prior (DCP) for removing haze in one
single image. It is observed that for non-sky outdoor areas, the intensity
value of the dark channels of a haze-free image is close to zero. Dark
channels are defined as the minimal value of three color channels of all
pixels within a local patch, which can be defined as follows:

$J^{dark}(x)=\min_{y\in\Omega(x)}\left(\min_{c\in\\{r,g,b\\}}J^{c}(y)\right),$
(3)

where $J^{c}$ is one color channel of $J$, $\\{r,g,b\\}$ is the three color
channels, $\Omega(x)$ is the local patch centered at pixel $x$. The statistics
of the haze-free outdoor images suggests that

$J^{dark}(x)\approx 0.$ (4)

DCP Atmospheric Light Estimation. In DCP, the global atmospheric light is
estimated from the most haze opaque areas, which are the pixels with the
highest dark channel value. More specifically, it ranks the dark channel
pixels by intensity values from high to low, then picks the top $0.1\%$
pixels. For these pixels, the one with the highest intensity value in the
original hazy image $I$ is used as the estimation of atmospheric light. As
shown in Fig. 1, the patch size will influence the final atmospheric light
estimation dramatically. A small patch size will make the process more subject
to noise. On the contrary, a large patch size will lead to inaccurate
estimation.

Fig. 1: Areas used to estimate the global atmospheric light, percentage=1%.
From left to right are the hazy image, areas using local patch size $3\times
3$ and local patch size $32\times 32$ respectively.

DCP Transmission Map Estimation. Despite the fairly good dehazing performance
of DCP, the same problem of patch size selection holds in the computation of
the transmission map. As shown in Fig. 2, details of hazy images can be
preserved in the refined transmission map when a small patch size is applied,
but it also introduces undesirable noise. When a large patch size is used, the
transmission map becomes smoother, while the sharpness is reduced. Therefore,
an appropriate patch size is also of great importance for a high-quality
estimation of the transmission map.

Fig. 2: Transmission map estimation. From left to right are the hazy image,
refined transmission map estimation using local patch size $3\times 3$ and
local patch size $32\times 32$ respectively.

###  2.2 Pyramid Fusion Dark Channel Prior (PF-DCP)

In Section 2.1, we introduce the basic ideas of DCP. The problems related to
the patch size are exposed yet not solved. To address this issue, we propose
the pyramid fusion dark channel prior (PF-DCP) based on the original DCP.

Fig. 3: The architecture of Pyramid Fusion Dark Channel Prior.

The key steps of PF-DCP are shown in Fig. 3. Note that we only show three
levels of the image pyramid, in practice, we generate more levels of low-
resolution images as long as their width and height are both above the patch
size. At each level, the hazy image performs steps similar to DCP, including
computing dark channels, estimating atmospheric light, estimating and refining
transmission maps. The above procedures are similar to those of DCP. Later,
after obtaining the transmission map at each level the lower-resolution
transmission map is firstly up-sampled to the same size as the higher-
resolution one, then we perform the weighted fusion among them. The process is
conducted iteratively until the highest resolution, i.e., the original image
size. The final transmission map will be further refined to recover the haze-
free image.

Image Pyramid Generation. We generate a multi-level image pyramid starting
from the original hazy image. We down-sample the image to half size until the
size is less than the patch size. The level of the pyramid depends on the
original size. In practice, it is normally within ten levels. After the
refined transmission map estimation is computed at each level, it is up-
sampled and restored iteratively to the original size. The transmission map
with the original size is then used to recover the haze-free image.

PF-DCP Atmospheric Light Estimation. To address the atmospheric light
estimation problem of the original DCP algorithm, we estimate it using dark
channels from all pyramid levels, before the computation of the transmission
maps. In this paper, we choose the maximal estimation value among all levels
as the final atmospheric light estimation.

PF-DCP Transmission Map Fusion. Due to the absence of the evaluation metrics
without ground truth, it is infeasible to perform automatic patch size
selection. Therefore, we apply a fusion-based algorithm to incorporate multi-
scale information. In this paper, we propose to perform an iterative weighted
fusion of the transmission maps at each level and dehaze the original image
via the final transmission map. In this way, the details of the original
images are kept as much as possible while the final transmission map preserves
multi-scale information. The reasons that we adopt a weighted fusion strategy
between the transmission map lie in three aspects:

  * •

The up-sampled transmission map from lower resolution is structurally similar
to the transmission map computed in higher resolution. Thus, performing a
fusion between them is practical and technically sound.

  * •

Different weight values offer much flexibility. If a smaller weight is set for
the transmission map from the lower resolution, the final transmission map
will contain more details of the scene. On the contrary, if a higher weight is
set, the final transmission map might lose some details, but noises are
reduced. Meanwhile, at least, part of the contribution of the final
transmission map comes from the highest resolution, so details are preserved
to some extent.

  * •

The weight fusion is computationally efficient and can be an alternative way
of different patch size selection. Specifically, a smaller patch size
maintains more details and a larger patch size reduces noise in the
transmission map. However, a larger patch size means a much slower computation
speed. On the contrary, changing the weight has nearly no impact on the
execution time, and has an even better performance compared to changing the
patch size.

Fused Transmission Map Refinement. To further improve the quality of the
results, we pay extra attention to the sharpness of the transmission map.
Inspired by DCP, we apply guided image filtering [18] to further refine the
fused transmission map. The guided filter assumes a local linear relationship
between the guiding image, which is the hazy image, and the output, which is
the refined transmission map. Utilizing the structural information of the hazy
image, the artifacts during the process of up-sampling and fusion are greatly
reduced.

Haze-free Image Recovery. After obtaining the final transmission maps of the
hazy images, we are able to easily recover the haze-free images using Eq. 2.
However, for the most haze opaque areas where the transmission $t$ is
extremely small, dividing $I(x)-A$ by $t$ is prone to be influenced by noise
[3]. Thus, we adopt a similar strategy to DCP to include a lower bound $t_{0}$
to make the algorithm more robust. The final formula to recover $J$ is

$J(x)=\frac{I(x)-A}{\max(t(x),t_{0})}+A.$ (5)

##  3 Experiments

###  3.1 Experimental Settings

By default, the patch size is $15\times 15$ for all pyramid levels. Top
$0.1\%$ haze opaque pixels are used for global atmospheric light estimation
and the lower bound of the transmission is $0.1$. The guided filter [18] is
adopted for transmission map refinement. To generate the image pyramid, we
discard half of the rows and columns directly. For the image up-sampling, we
use the nearest-neighbor interpolation. In our in-house experiments, there is
no obvious performance difference by using other interpolation algorithms,
such as bilinear, bicubic and nearest-neighbor. Using the extra training set
of RESIDE [19], the best fusion weight of two transmission maps for every
consecutive pairs, the lower and higher resolution level, is $4:1$ for indoor
scenes $80:1$ for outdoor scenes. Compared to the original DCP, the
computation cost of PF-DCP increases by $30\%$ on average.

For the parameter selection, we use images from Indoor Training Set (ITS) and
Outdoor Training Set (OTS) of RESIDE [19]. For testing, we evaluate the
performance of the algorithm on the Synthetic Objective Testing Set (SOTS) of
RESIDE, which contains 500 labeled indoor hazy images and 500 outdoor ones,
and OHaze [20] for real hazy images. We use two commonly-used metrics for
evaluation in single image dehazing: Peak Signal-to-Noise Ratio (PSNR) and
Structural Similarity (SSIM).

###  3.2 Experimental Results

Fusion Weight. Since we have already generated an image pyramid and calculated
transmission map at each level in PF-DCP, setting a higher weight value at the
desired level of the image pyramid shall achieve a similar effect of selecting
the desired patch size. Through experiments, we find out that fixating the
patch size and changing the fusion weight achieves better performance both in
speed and in quality than fixating fusion weight and changing patch size.

For simplification, we fix the weight value between every two levels in the
image pyramid. Fig. 4 shows the results of using different fusion weights. The
performance of outdoor scenes is much better than that of indoor ones. It is
not surprising because DCP is for outdoor scenes in theory. For indoor scenes,
they are sensitive to the change of fusion weights and achieves the best
performance at a relatively small weight value. For outdoor scenes, they are
the opposite. The plausible reasons are that there are more details and the
relative depth change is large in a local patch of the indoor scenes.

Fig. 4: Results of different fusion weights on RESIDE ITS and OTS.

Comparison with state-of-the-art methods. The proposed PF-DCP is compared
against traditional prior-based approaches, DCP [3], BCCR [4], CAP [5], and
NLD [6]. [19].

(a)

(b)

(c)

(d)

(e)

(f)

(g)

Fig. 5: Dehazed images of OHaze (first 2 rows) and RESIDE SOTS (remaining 4
rows). (a) Hazy images. The restored images using (b) BCCR (c) CAP (d) NLD (e)
DCP (f) PF-DCP with fusion weight $4:1$(indoor) and $80:1$(outdoor) (g) Ground
truth.

As we could see from Fig. 5, BCCR introduces large color distortion, and the
same problem holds for NLD. The results are also reflected in the quantitative
analysis. As shown in Table 1, the PSNR scores for both BCCR and NLD are low
due to the inauthentic color. For CAP, the colors are close to those in
reality but there is still a lot of haze remaining after dehazing. Thus, CAP
achieves a fairly good yet not satisfactory performance in PSNR and SSIM. As
for the original DCP, there are apparent halo artifacts around the edges of
objects, which comes from the inaccurate transmission rate estimation.
Besides, the dehazed images are dark because of the faulty atmospheric light
estimation. On the contrary, utilizing the multi-scale information, PF-DCP
maintains the authentic color and has little halo artifacts. Therefore, PSNR
and SSIM scores for PF-DCP get largely improved compared to the original DCP.

BCCR [4] | CAP [5] | NLD [6] | AOD [9] | DCPDN [15] | GFN [10] | EPDN [14] | DCP [3] | PF-DCP
---|---|---|---|---|---|---|---|---
16.30/0.80 | 20.72/0.87 | 17.68/0.81 | 19.78/0.87 | 17.98/0.86 | 21.81/0.86 | 22.68/0.90 | 17.82/0.86 | 23.07/0.91

Table 1: Average PSNR and SSIM on the RESIDE SOTS.

In addition to traditional methods, we also perform a quantitative evaluation
between PF-DCP with contemporary deep learning methods on RESIDE SOTS. As
shown in Table 1, PF-DCP manages to outperform listed popular deep learning
methods in recent years. It is true that the proposed PF-DCP is still inferior
to the most recent deep learning methods [21], however, PF-DCP has its own
advantages. First, PF-DCP provides reasonably good performance. Second, it
does not require computationally intensive model training. Third, it is
trackable and explainable rather than a black box. We believe the advantages
of PF-DCP reveal great potentials for single image dehazing in real
applications.

##  4 Conclusion

In this paper, we have presented an effective algorithm, pyramid fusion dark
channel prior (PF-DCP) by introducing multi-scale information of hazy images
based on the original DCP. Through experiments, we show that fusing the
transmission maps of different image scales can boost the performance of haze
removal greatly. Meanwhile, the weighted fusion strategy provides users an
alternative to the patch size selection, which is more robust and
computationally efficient.

## References

  * [1] Nicolas Hautiére, Jean-Philippe Tarel, Jean Lavenant, and Didier Aubert,  “Automatic fog detection and estimation of visibility distance through use of an onboard camera,”  Machine Vision and Applications, vol. 17, no. 1, pp. 8–20, 2006\.
  * [2] Robby T Tan,  “Visibility in bad weather from a single image,”  in 2008 IEEE Conference on Computer Vision and Pattern Recognition. IEEE, 2008, pp. 1–8.
  * [3] Kaiming He, Jian Sun, and Xiaoou Tang,  “Single image haze removal using dark channel prior,”  IEEE transactions on pattern analysis and machine intelligence, vol. 33, no. 12, pp. 2341–2353, 2010.
  * [4] Gaofeng Meng, Ying Wang, Jiangyong Duan, Shiming Xiang, and Chunhong Pan,  “Efficient image dehazing with boundary constraint and contextual regularization,”  in Proceedings of the IEEE international conference on computer vision, 2013, pp. 617–624.
  * [5] Qingsong Zhu, Jiaming Mai, and Ling Shao,  “A fast single image haze removal algorithm using color attenuation prior,”  IEEE transactions on image processing, vol. 24, no. 11, pp. 3522–3533, 2015.
  * [6] Dana Berman, Shai Avidan, et al.,  “Non-local image dehazing,”  in Proceedings of the IEEE conference on computer vision and pattern recognition, 2016, pp. 1674–1682.
  * [7] Bolun Cai, Xiangmin Xu, Kui Jia, Chunmei Qing, and Dacheng Tao,  “Dehazenet: An end-to-end system for single image haze removal,”  IEEE Transactions on Image Processing, vol. 25, no. 11, pp. 5187–5198, 2016.
  * [8] Wenqi Ren, Si Liu, Hua Zhang, Jinshan Pan, Xiaochun Cao, and Ming-Hsuan Yang,  “Single image dehazing via multi-scale convolutional neural networks,”  in European conference on computer vision. Springer, 2016, pp. 154–169.
  * [9] Boyi Li, Xiulian Peng, Zhangyang Wang, Jizheng Xu, and Dan Feng,  “Aod-net: All-in-one dehazing network,”  in Proceedings of the IEEE international conference on computer vision, 2017, pp. 4770–4778.
  * [10] Wenqi Ren, Lin Ma, Jiawei Zhang, Jinshan Pan, Xiaochun Cao, Wei Liu, and Ming-Hsuan Yang,  “Gated fusion network for single image dehazing,”  in Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, 2018, pp. 3253–3261.
  * [11] Gary Marcus,  “Deep learning: A critical appraisal,”  arXiv preprint arXiv:1801.00631, 2018.
  * [12] Haoran Xu, Jianming Guo, Qing Liu, and Lingli Ye,  “Fast image dehazing using improved dark channel prior,”  in 2012 IEEE International Conference on Information Science and Technology. IEEE, 2012, pp. 663–667.
  * [13] Jinshan Pan, Deqing Sun, Hanspeter Pfister, and Ming-Hsuan Yang,  “Blind image deblurring using dark channel prior,”  in Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, 2016, pp. 1628–1636.
  * [14] Yanyun Qu, Yizi Chen, Jingying Huang, and Yuan Xie,  “Enhanced pix2pix dehazing network,”  in Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, 2019, pp. 8160–8168.
  * [15] He Zhang and Vishal M Patel,  “Densely connected pyramid dehazing network,”  in Proceedings of the IEEE conference on computer vision and pattern recognition, 2018, pp. 3194–3203.
  * [16] Yu Dong, Yihao Liu, He Zhang, Shifeng Chen, and Yu Qiao,  “Fd-gan: Generative adversarial networks with fusion-discriminator for single image dehazing.,”  in AAAI, 2020, pp. 10729–10736.
  * [17] Sungmin Lee, Seokmin Yun, Ju-Hun Nam, Chee Sun Won, and Seung-Won Jung,  “A review on dark channel prior based image dehazing algorithms,”  EURASIP Journal on Image and Video Processing, vol. 2016, no. 1, pp. 4, 2016.
  * [18] Kaiming He, Jian Sun, and Xiaoou Tang,  “Guided image filtering,”  in European conference on computer vision. Springer, 2010, pp. 1–14.
  * [19] Boyi Li, Wenqi Ren, Dengpan Fu, Dacheng Tao, Dan Feng, Wenjun Zeng, and Zhangyang Wang,  “Reside: A benchmark for single image dehazing,”  arXiv preprint arXiv:1712.04143, vol. 1, 2017.
  * [20] Codruta O Ancuti, Cosmin Ancuti, Radu Timofte, and Christophe De Vleeschouwer,  “O-haze: a dehazing benchmark with real hazy and haze-free outdoor images,”  in Proceedings of the IEEE conference on computer vision and pattern recognition workshops, 2018, pp. 754–762.
  * [21] Yuanjie Shao, Lerenhan Li, Wenqi Ren, Changxin Gao, and Nong Sang,  “Domain adaptation for image dehazing,”  in Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, 2020, pp. 2808–2817.