# Temporally Modulated One-Dimensional Leaky-Wave Holograms

Amrollah Amini  Iran University of Science and Technology, School of
Electrical Engineering, Tehran, 1684613114, Iran  Homayoon Oraizi
<EMAIL_ADDRESS>
###### Abstract

Spatio-temporally modulated impedance surfaces can be good candidates for
generation of radiating waves with arbitrary eigenstates by breaking momentum
and energy conservations. Here, we present a theoretical framework based on
the holographic technique and generalized Floquet-wave expansion to analyze
spatio-temporally modulated impedance surfaces. The holographic technique
estimates the required impedance distribution to achieve the desired momentum.
Injecting temporal modulation deviates the eigenvalues and changes the
radiation frequency. Using the proposed analytical model, the eigenvalues can
be calculated accurately in the presence of space and time modulations.
Consequently, it is possible to predict the propagation mechanism of bounded
and radiation states. It has been shown that, imposition of temporal
modulation causes the Doppler-shift effect and nonreciprocal responses in the
hologram. By plotting the antenna dispersion diagram, and observing the
asymmetric displacement of dispersion curve due to temporal modulation, the
system nonreciprocity can be verified. The beam scanning properties of these
structures have also been investigated.

## Introduction

Recently, planar nonreciprocal devices without the use of bulky ferromagnetic
materials have received much attention [1, 2] due to their subwavelength
thicknesses and low cost manufacturing processes. A suitable solution for
implementing such devices is to exploit spatio-temporally modulated
structures. This method can cover a wide range of frequencies from microwaves
[3] to optical [4, 5] regimes. The use of spatio-temporally modulated
structures is not limited to electromagnetic devices, but can also be applied
to manipulate acoustic waves [6]. A wide variety of potential applications can
be envisaged for space-time structures, such as isolators [7], frequency
translators [8], circulators [9], and nonreciprocal phase shifters [10]. In
the context of antenna engineering several works have also been proposed to
realize nonreciprocal radiators [11, 12, 13]. A mixer-duplexer leaky-wave
antenna is presented in the literature using temporally modulated microstrip
lines [12]. An array of varactors incorporated in transmission lines was used
to inject temporal modulation. A silicon-compatible leaky-wave antenna has
been designed based on graphene impedance surfaces [13]. Nonreciprocity and
Doppler-shift effects are obtained by applying temporal modulation through
gating pads located adjacent to graphene surfaces. These structures can be
used for bio-sensing, imaging, and inter-chip communication applications [13].

An appropriate approach to implement spatio-temporally modulated structures is
to use dynamic metasurfaces, which have several advantages such as ultra-thin
thicknesses, low cost and high flexibility. Metasurfaces are mainly exploited
in transmissive [14, 15, 16, 17, 18, 19], reflective [20, 21, 22, 23], or
leaky-wave [24, 25] modes. The main drawback of microwave metareflectors and
transmitarrays is their protruded feeds which make them unsuitable for low-
profile integrated systems. This shortcoming can be resolved by using leaky-
wave metasurfaces which are generally fed by embedded monopole structures.
Owing to the above-mentioned attractive features, various methods have been
recently proposed for the design and analysis of leaky-wave metasurfaces. The
holographic technique [26, 27, 28], aperture field estimation (AFE) method
[29, 30, 31], Floquet-wave expansion model [32, 33] and the method of moments
(MoM) [34, 35] are among the proposed methods of designing metasurface-based
leaky-wave antennas. The holographic and AFE techniques can be considered as
synthesis methods for the estimation of arrangement of meta-atoms to achieve
the desired wavefronts, whereas the MoM and Floquet-wave expansion are
analysis frameworks for the determination of exact solutions of leaky-wave
modes generated by the synthesized metasurfaces. Therefore, for a
comprehensive study of leaky-wave metasurfaces, we need to combine the
analysis and synthesis methods.

In this paper, the combination of holographic technique and generalized
Floquet-wave expansion is proposed to analyze the temporally modulated one-
dimensional leaky-wave metasurfaces. The scalar impedance boundary condition
is used to model the proposed hologram. The generalized Floquet-wave expansion
was first proposed by Cassedy [36] to attain the mode characteristics of
temporally modulated surface waveguides. In our work, this method is used for
impedance surfaces to explain the radiation mode and the leakage mechanism in
the presence of space-time modulation. As an example, using the holographic
technique, the surface impedance is designed to radiate in the direction of
$30^{\circ}$ at 18 GHz. The effect of temporal injection is studied and the
radiation characteristics of antenna in the transmission and reception modes
have been investigated.

## Spatio-temporally modulated impedance surfaces

A convenient way to implement holographic metasurfaces in the microwave regime
is to modulate impedance surfaces. For an impedance surface located at $z=0$
(see Fig. 1), the tangential components of electric ($\vec{E}_{t}$) and
magnetic ($\vec{H}_{t}$) fields are related through the impedance boundary
condition as follows [37]:

$\vec{E}_{t}(x,y)=\underline{\underline{Z}}_{s}(x,y)\cdot\hat{z}\times\vec{H}_{t}(x,y)$
(1)

where $\underline{\underline{Z}}_{s}$ is the tensorial impedance indicating
anisotropic boundary conditions and $\hat{z}$ denotes the unit vector along
the z-axis. For scalar (isotropic) impedance surfaces, the boundary condition
is simplified as:

$\vec{E}_{t}(x,y)=Z_{s}(x,y)\hat{z}\times\vec{H}_{t}(x,y)$ (2)

Figure 1 shows the conceptual representation of an impedance boundary
condition. For spatially modulated surfaces, dielectric slabs covered by
periodic subwavelength patches can be used for the realization of impedance
boundary conditions. Generally, impedance surfaces support both TM and TE
surface modes. Note that for the TM (TE) mode the surface wave has transverse
magnetic (electric) field components in both propagation and $\hat{z}$
directions. In [38] it has been shown that if the surface impedance is
modulated sinusoidally in the propagation direction, the power will leak from
the surface at a certain angle which is proportional to the impedance
periodicity. In this case, the surface impedance acts as a leaky-wave antenna.
However, this structure has several characteristics such as linearity and
reciprocity. These properties may be rectified by adding temporal variations
to the impedance surface. In spatio-temporally modulated surfaces, the
impedance boundary condition changes simultaneously in space and time domains.
In the general case, the impedance boundary condition can be expressed as
follows:

$Z_{s}(x,y;t)=\sum_{n}Z_{n}e^{-jnKs(x,y;t)}$ (3)

where $Ks$ is the modulation phase and represents the local periodic
variations in the impedance boundary condition. The coefficient $Z_{n}$ can be
a complex number and determines the variation depth of surface impedance. For
a monochromatic modulation, equation (3) is simplified as:

$Z_{s}(x,y;t)=Z_{0}+Z_{-1}e^{jKs(x,y;t)}+Z_{+1}e^{-jKs(x,y;t)}$ (4)

In this case we may define the modulation frequencies in space and time
domains, which are indeed the spatial and temporal derivatives of $Ks(x,y;t)$,
respectively:

$\vec{\beta}_{p}=\nabla_{x,y}Ks(x,y;t)=\frac{\partial}{\partial
x}Ks(x,y;t)\hat{x}+\frac{\partial}{\partial y}Ks(x,y;t)\hat{y}$ (5)
$\omega_{p}=-\frac{\partial}{\partial t}Ks(x,y;t)$ (6)

The frequencies $\vec{\beta}_{p}$ and $\omega_{p}$ are actually space and time
domain pumping frequencies, respectively, that cause power to be coupled from
surface mode to higher order Floquet modes. If the pumping frequencies are
equal to zero, the power coupling will not take place and we will have a pure
surface mode. Note that, $Ks$ can be an arbitrary function of (x, y) in the
cartesian coordinate system. In this paper, we assume that $Ks$ varies only
along the $\hat{x}$ direction and is constant along the $\hat{y}$ axis. Figure
2 shows an example of spatio-temporally modulated impedance surface in
different time slots for $Ks(x;t)=\beta_{p}x-\omega_{p}t$ and $Z_{-1}=Z_{+1}$.
In the next section, we will see that the impedance synthesized by the
holographic method will be similar to that presented in Fig.2. The progressive
wave-like surface pattern travels with speed of
$\nu_{p}=-\partial_{t}Ks(x)/(\nabla_{x}Ks(x).\hat{x})=\omega_{p}/\beta_{p}$
along the $\hat{x}$ axis. Note that, since the boundary condition of the
problem varies with time, we must use the wave equation in the time domain to
analyze such structures. Thus, for the half-space above the surface the wave
equation can be written as:

$(\frac{\partial^{2}}{\partial x^{2}}+\frac{\partial^{2}}{\partial
z^{2}}-\frac{1}{c^{2}}\frac{\partial^{2}}{\partial t^{2}})\psi(x,z;t)=0$ (7)

where c is the speed of light in vacuum. Note that, for a single frequency
system without time modulation, assuming that the wave function depends on
$e^{j\omega t}$, the wave equation results in the Helmholtz equation. However,
in the presence of a time-varying boundary condition and the effects of higher
order Floquet harmonics, the Helmholtz equation is not valid. The periodicity
of boundary condition in both the time and space domains requires that the
generalized Floquet-wave theory [36] be applied. In this case, the wave
function ($\psi(x,z;t)$) can be written in terms of higher order modes [36]:

$\psi(x,z;t)=\sum_{n}\psi_{n}e^{j(\omega+n\omega_{p})t}e^{-j(\kappa+n\beta_{p})x}e^{-jk_{zn}z}$
(8)

where $\omega$ and $\kappa$ are the fundamental mode frequencies in time and
space domains, respectively. The coefficient $\psi_{n}$ indicates the harmonic
amplitude and $jk_{zn}$ is the propagation constant of n-indexed harmonic in
the $\hat{z}$ direction. Using (7), $k_{zn}$ can be determined as follows:

$k_{zn}=\sqrt{\frac{(\omega+n\omega_{p})^{2}}{c^{2}}-(\kappa+n\beta_{p})^{2}}$
(9)

In equation (8) $\omega_{p}$ and $\beta_{p}$ can be considered as degrees of
freedom, of which the values may be determined by the designer. To estimate
$\beta_{p}$, we can use holographic theory for spatially modulated surfaces.

## Holographic technique for leaky-wave metasurfaces

An effective approach for the design and synthesis of leaky-wave modulated
metasurface antennas is the holographic technique, which was first proposed by
Sievenpiper et al. [26, 27]. In this technique, for a hologram, the required
distribution of surface impedance to generate an object wave with desired
direction (called $\psi_{obj}$) can be obtained as:

$Z_{s}(x,y)=jX_{0}[1+M\times\Re\\{\psi_{ref}^{*}(x,y)\psi_{obj}(x,y)\\}]$ (10)

where $\psi_{ref}$ is the reference wave excited by the electromagnetic source
and asterisk ($*$) denotes the complex conjugate operation. This gradient
distribution of surface impedance plays the role of holographic interferogram
building up the object wave as the reference wave excites it. In (10), $X_{0}$
and M are average surface reactance and modulation depth, respectively. For
leaky-wave holograms the modulated depth must be selected small enough (M <
0.6) in order to achieve the specified directive radiations [39]. If the
radiated wave is supposed to be directed in the $\theta=\theta_{0}$ zenith
angle, the corresponding wave (object wave) can be represented by:

$\psi_{obj}(x)=e^{-jk_{d}\sin\theta_{0}x}$ (11)

where $k_{d}$ indicates the wave-number at the design frequency $f=f_{d}$. For
a planar wave-front propagating along the $\hat{x}$ axis, the reference wave
can be defined as [39]:

$\psi_{ref}(x)=e^{-jk_{d}\sqrt{1+(\frac{X_{0}}{\eta_{0}})^{2}}x}$ (12)

where $\eta_{0}=120\pi$ is the impedance of free space. Substituting (11) and
(12) in (10) yields

$Z_{s}(x)=jX_{0}[1+M\times\cos(k_{d}(\sqrt{1+(\frac{X_{0}}{\eta_{0}})^{2}}-\sin\theta_{0})x)]$
(13)

Considering the modulation phase and comparing it with the impedance
distribution in the previous section, we can conclude that:

$Ks(x)=k_{d}(\sqrt{1+(\frac{X_{0}}{\eta_{0}})^{2}}-\sin\theta_{0})x$ (14)
$Z_{0}=jX_{0},\quad Z_{+1}=Z_{-1}=j\frac{X_{0}M}{2}$ (15)
$\beta_{p}=k_{d}(\sqrt{1+(\frac{X_{0}}{\eta_{0}})^{2}}-\sin\theta_{0})$ (16)

Time dependence can be applied to the modulation phase in different ways.
However, for simplicity, we assume that the temporal variation of the
modulation is linear. Therefore, for dynamic structures, $Ks(x;t)$ may be
redefined as follows:

$Ks(x;t)=k_{d}(\sqrt{1+(\frac{X_{0}}{\eta_{0}})^{2}}-\sin\theta_{0})x-\omega_{p}t$
(17)

## Calculation of propagation constant for asymptotic case

If the surface wave excited by the launcher is in the form of $TM_{0}$ mode,
the x and z components of electric field ($E_{x}$ and $E_{z}$) and the y
component of magnetic field ($H_{y}$) may contribute to the propagation along
the impedance boundary condition. Using the generalized Floquet-wave expansion
we have:

$E_{x}(x,z;t)=\sum_{n}E_{n}e^{j(\omega+n\omega_{p})t}e^{-j(\kappa+n\beta_{p})x}e^{-jk_{zn}z}$
(18)
$H_{y}(x,z;t)=\sum_{n}H_{n}e^{j(\omega+n\omega_{p})t}e^{-j(\kappa+n\beta_{p})x}e^{-jk_{zn}z}$
(19)

From Maxwell’s equations the z component of electric field can be obtained:

$E_{z}(x,z;t)=-\sum_{n}\frac{(\kappa+n\beta_{p})H_{n}}{(\omega+n\omega_{p})\epsilon_{0}}e^{j(\omega+n\omega_{p})t}e^{-j(\kappa+n\beta_{p})x}e^{-jk_{zn}z}$
(20)

In the above expansions, $\kappa$, $E_{n}$ and $H_{n}$ are unknown parameters
and must be determined. The n-indexed Floquet mode frequencies may also be
defined as follows:

$\beta^{(n)}=\kappa+n\beta_{p}$ (21) $\omega^{(n)}=\omega+n\omega_{p}$ (22)

The Eigenmode analysis is an appropriate method to calculate the unknowns
($\kappa$, $E_{n}$ and $H_{n}$). Substituting (18) and (19) in the impedance
boundary condition, and comparing the terms with identical modulation phases,
leads to

$E_{n}=-jX_{0}[H_{n}+\frac{M}{2}(H_{n+1}+H_{n-1})]$ (23)

Furthermore, the Maxwell equation imposes the following restriction:

$\frac{E_{n}}{H_{n}}=\frac{k_{zn}}{(\omega+n\omega_{p})\epsilon_{0}}$ (24)

Combining (23) and (24) yields:

$H_{n+1}+D_{n}H_{n}+H_{n-1}=0\quad n=0,\pm 1,\pm 2,...$ (25)

where:

$D_{n}=\frac{2}{M}[1-\frac{jk_{zn}}{X_{0}(\omega+n\omega_{p})\epsilon_{0}}]$
(26)

Equation (25) relates the magnitude of n-indexed harmonic to (n+1) and
(n-1)-indexed modes, which form a set of infinite equations with infinite
unknowns. These equations resemble the dispersion relations in the pure-space
modulated case [38], except that $\omega^{(n)}=\omega+n\omega_{p}$ is used
instead of $\omega$. In order to solve this system of equations, we must
truncate the number of modes for analysis.

For the asymptotic case, if the modulation index tends to zero ($M\rightarrow
0$), the numerator of fraction in equation (26) needs to be zero to achieve
nontrivial solutions. In this case we have:

$1-\frac{j\sqrt{\frac{(\omega+n\omega_{p})^{2}}{c^{2}}-(\kappa+n\beta_{p})^{2}}}{X_{0}(\omega+n\omega_{p})\epsilon_{0}}=0$
(27)

Therefore, the set of asymptotic curves for $\kappa$ can be obtained from the
following equation:

$\kappa=-nk_{d}(\sqrt{1+(\frac{X_{0}}{\eta_{0}})^{2}}-\sin\theta_{0})\pm\frac{\omega+n\omega_{p}}{c}\sqrt{1+(\frac{X_{0}}{\eta_{0}})^{2}}\quad
n=0,\pm 1,\pm 2,...$ (28)

Figure 3 shows the asymptotic dispersion curves for the cases of $f_{p}=0$,
$f_{p}=1GHz$ and $f_{p}=2GHz$. The solid and dashed curves represent forward
and backward solutions, respectively. These curves can be used as initial
guesses for the calculation of eigenstates of the modulated case. Observe
that, injection of temporal variations into the modulated impedance function
imposes some asymmetry on the dispersion curves, leading to a nonreciprocal
response. As shown in Fig. 3, in the presence of temporal modulation, the
Brillouin diagram is parallel to the line with slope of $\nu_{p}/c$ and the
period of Brillouin zone (the distance between adjacent red points in Fig. 3)
is $\sqrt{\beta_{p}^{2}+(\omega_{p}/c)^{2}}$. Observe that increasing the
temporal pumping frequency results in the increase of the degree of asymmetry.
This property will be described in more detail in the following section.

## Determination of propagation constant in the presence of space-time
modulation

As discussed in the previous section, The infinite system of equations in (25)
should be solved for the determination of the exact values of eigenstates in
the presence of space-time modulation, which is achieved by safely truncating
the number of modes without generating of a significant error. In reference
[40] it has been shown that considering only three Floquet modes (namely
$n=0,\pm 1$) results in sufficiently accurate solutions. Therefore, the
following equations will be formed:

$H_{0}+D_{-1}H_{-1}=0$ (29) $H_{1}+D_{0}H_{0}+H_{-1}=0$ (30)
$D_{1}H_{1}+H_{0}=0$ (31)

Combining (29)-(31), yields:

$(D_{0}-\frac{1}{D_{-1}}-\frac{1}{D_{1}})H_{0}=0$ (32)

To obtain nontrivial solutions, it is required to have:

$D_{0}-\frac{1}{D_{-1}}-\frac{1}{D_{1}}=0$ (33)

Therefore, we can rewrite the dispersion equation as:

$\frac{4}{M^{2}}[1-\frac{j\sqrt{(\omega/c)^{2}-\kappa^{2}}}{X_{0}\omega\epsilon_{0}}]=\frac{1}{1-\frac{j\sqrt{(\omega-\omega_{p})^{2}/c^{2}-(\kappa-\beta_{p})^{2}}}{X_{0}(\omega-\omega_{p})\epsilon_{0}}}+\frac{1}{1-\frac{j\sqrt{(\omega+\omega_{p})^{2}/c^{2}-(\kappa+\beta_{p})^{2}}}{X_{0}(\omega+\omega_{p})\epsilon_{0}}}$
(34)

Its solution determines the complex values of wave-number at a given
frequency. Figure 4a shows the dispersion diagram of hologram for the case of
$f_{p}=0\,GHz$. The design frequency is 18 GHz, and the beam angle at 18 GHz
is set at $\theta_{0}=30^{\circ}$. The average impedance ($X_{0}$) and
modulation index (M) are selected as $0.8\eta_{0}$ and 0.4. Therefore, the
corresponding impedance may be obtained as:

$Z_{s}(x)=j0.8\eta_{0}[1+0.4\cos(377.25(\sqrt{1+0.8^{2}}-0.5)x)]=j301.59[1+0.4\cos(294.49x)]$
(35)

The total magnetic field at 18 GHz is also plotted in Fig. 4b, which show
that, in the transmission mode, the coupling occurs from the surface to
radiation mode with the radiation angle of 30 degrees. Also, in the reception
mode, the maximum coupling occurs at the incoming angle of 30 degrees, which
confirms the accuracy of the design. In holographic antennas, the coupling
strength depends on the modulation depth (M). Figure 5 shows the total fields
for different values of modulation depths. Observe that for lower values of
$M$, most of the power is bound to the surface. The power coupling is
increased by enhancing the modulation depth.

Injecting temporal modulation into the designed hologram, imposes some
modifications to the dispersion diagram. Figure 6 shows the dispersion curves
for the case of $f_{p}=2\,GHz$. Its impedance distribution has the following
form:

$Z_{s}(x;t)=j301.59[1+0.4\cos(294.49x-4\pi\times 10^{9}t)]$ (36)

Observe that, in the presence of temporal modulation, the curves corresponding
to negative and positive harmonics move down and up respectively. This
displacement is proportional to the pumping frequency ($f_{p}$). As the
pumping frequency increases, the dispersion curve for $n=-1$ tends towards the
light-line, indicating that the radiation mode tends towards the end-fire
direction. Also for $n=+1$, the radiation curve gets closer to the vertical
axis. To investigate nonreciprocity, we first consider the antenna operation
in transmission mode. Suppose that the excitation frequency of the antenna is
18 GHz, which operates at the fundamental mode (surface wave). According to
the dispersion diagram in Fig. 6a, its corresponding phase constant is equal
to $\beta^{(0)}=488\,(rad/m)$. Furthermore, using (21) the phase constant for
radiation mode ((-1)-indexed mode) can be obtained as
$\beta^{(-1)}=\beta^{(0)}-\beta_{p}\approx 185\,(rad/m)$. Observe in the
dispersion diagram that the corresponding frequency at this phase constant is
16 GHz (namely $f^{(-1)}=f^{(0)}-f_{p}$), which shows the Doppler-shift effect
in the proposed antenna. According to the calculated phase constant, the beam
orientation is directed at 35 degrees. In Fig. 7 the analytical results of
total magnetic field and Floquet-modes of $n=0,\pm 1$ are plotted separately.
Figure 7 indicates that a strong coupling can be observed between 0-indexed
and (-1)-indexed Floquet modes. In this case, most of the power is divided
into the fundamental and (-1)-indexed modes, but the (+1)-indexed mode has a
small portion of the power. The radiation mode is oriented in the direction of
35 degrees. The blue and red arrows in Fig. 7a indicate the beam directions
for pure-space and space-time modulated cases, respectively. When the antenna
is in the reception mode, for the 16 GHz incoming wave, the corresponding
phase constant is entirely different from that of the transmission case. Using
(22), for (+1)-indexed mode with frequency of 16 GHz, the 0-indexed mode
propagates at $f^{(0)}=f^{(+1)}-f_{p}=14GHz$. In this case, the phase constant
of (+1)-indexed mode is equal to $\beta^{(+1)}\approx-79\,(rad/m)$ and the
angle at which the maximum power is coupled to the surface is 14 degrees (as
shown in Fig. 8). Note that the surface is excited at 14 GHz instead of 16
GHz.

To summarize our analysis we can conclude that, for the antenna in the
transmission mode, the maximum coupling between the surface and radiation
modes occurs at the angle of 35 degrees, whereas for the reception case the
radiation and surface modes have the maximum coupling at the incoming angle of
14 degrees . Therefore, the antenna does not behave the same in the
transmission and reception modes, which confirms its nonreciprocal behavior.

## Beam scanning properties

Leaky-wave holograms can be suitable choices for automotive sensors, coherent
tomography, real-time spectrum analysis and tracking applications owing to
their beam scanning capabilities, which can be controlled by source frequency.
Referring to the dispersion diagram in Fig. 4, for a pure-space modulated
hologram, the corresponding phase constant increases by enhancing excitation
frequency. This means that, the radiation beam tends towards the end-fire
direction at higher frequencies. In temporally modulated holograms, beam
scanning can be achieved through another mechanism. In this case, the beam
direction varies by changing the pumping frequency ($f_{p}$) without the need
to change the source frequency. The calculated dispersion diagrams for
different values of $f_{p}$ are plotted in Fig. 9. Observe that for a given
phase constant, the (-1)-indexed curve tends towards the light-line as the
pumping frequency increases. The dispersion curve for (+1)-indexed mode moves
in opposite direction indicating that the incoming wave angle for reception
case tends towards the broadside direction. Figure 10 shows the calculated
fields for different values of $f_{p}$. Results show that the beam direction
changes from $\theta=30^{\circ}$ to $\theta=41^{\circ}$ with pumping frequency
increasing from 0 GHz to 4 GHz. It should be noted that the source frequency
is fixed at 18 GHz.

## Conclusion

In this paper, the generalized Floquet-wave expansion method is utilized to
accurately calculate the propagation characteristics of temporally modulated
holograms. The combination of holographic technique as a synthesis method and
Floquet-wave expansion can form a fully-analytical model for implementation of
leaky-wave metasurfaces. It has been shown that temporally modulated holograms
can be effectively used as nonreciprocal antennas. This nonreciprocity can be
obtained by properly displacing the dispersion curve so that the beam is
directed at the desired spherical angle. Owing to their embedded feeding
networks, leaky-wave antennas can be appropriate alternatives to nonreciprocal
metareflectors and transmitarrays as radiators of integrated transceivers.

## Author contributions statement

A. A. and H. O. developed the idea and theory, and wrote the manuscript. A. A.
did the simulations. H. O. analyzed the results and supervised the project.

## Additional information

Competing Interests: The authors declare no competing interests.

## References

  * [1] Estep, N. A., Sounas, D. L., Soric, J. & Alù, A.  Magnetic-free non-reciprocity and isolation based on parametrically modulated coupled-resonator loops.  _Nature Physics_ 10, 923–927, DOI: 10.1038/nphys3134 (2014).
  * [2] Sounas, D. L. & Alù, A.  Non-reciprocal photonics based on time modulation.  _Nature Photonics_ 11, 774–783, DOI: 10.1038/s41566-017-0051-x (2017).
  * [3] Taravati, S. & Kishk, A. A.  Space-time modulation: Principles and applications.  _IEEE Microwave Magazine_ 21, 30–56, DOI: 10.1109/MMM.2019.2963606 (2020).
  * [4] Shaltout, A. M., Shalaev, V. M. & Brongersma, M. L.  Spatiotemporal light control with active metasurfaces.  _Science_ 364, eaat3100, DOI: 10.1126/science.aat3100 (2019).
  * [5] Shaltout, A. M. _et al._ Spatiotemporal light control with frequency-gradient metasurfaces.  _Science_ 365, 374–377, DOI: 10.1126/science.aax2357 (2019).
  * [6] Fakheri, M. H., Rajabalipanah, H. & Abdolali, A.  Spatiotemporal binary acoustic metasurfaces.  _Phys. Rev. Applied_ 16, 024062, DOI: 10.1103/PhysRevApplied.16.024062 (2021).
  * [7] Yu, Z. & Fan, S.  Complete optical isolation created by indirect interband photonic transitions.  _Nature Photonics_ 3, 91–94, DOI: 10.1038/nphoton.2008.273 (2009).
  * [8] Wu, Z. & Grbic, A.  Serrodyne frequency translation using time-modulated metasurfaces.  _IEEE Transactions on Antennas and Propagation_ 68, 1599–1606, DOI: 10.1109/TAP.2019.2943712 (2020).
  * [9] Sounas, D. L., Caloz, C. & Alù, A.  Giant non-reciprocity at the subwavelength scale using angular momentum-biased metamaterials.  _Nature Communications_ 4, 2407, DOI: 10.1038/ncomms3407 (2013).
  * [10] Wang, X., Díaz-Rubio, A., Li, H., Tretyakov, S. A. & Alù, A.  Theory and design of multifunctional space-time metasurfaces.  _Phys. Rev. Applied_ 13, 044040, DOI: 10.1103/PhysRevApplied.13.044040 (2020).
  * [11] Taravati, S. & Eleftheriades, G. V.  Space-time medium functions as a perfect antenna-mixer-amplifier transceiver.  _Phys. Rev. Applied_ 14, 054017, DOI: 10.1103/PhysRevApplied.14.054017 (2020).
  * [12] Taravati, S. & Caloz, C.  Mixer-duplexer-antenna leaky-wave system based on periodic space-time modulation.  _IEEE Transactions on Antennas and Propagation_ 65, 442–452, DOI: 10.1109/TAP.2016.2632735 (2017).
  * [13] Correas-Serrano, D. _et al._ Nonreciprocal graphene devices and antennas based on spatiotemporal modulation.  _IEEE Antennas and Wireless Propagation Letters_ 15, 1529–1532, DOI: 10.1109/LAWP.2015.2510818 (2016).
  * [14] Khorasaninejad, M. _et al._ Metalenses at visible wavelengths: Diffraction-limited focusing and subwavelength resolution imaging.  _Science_ 352, 1190–1194, DOI: 10.1126/science.aaf6644 (2016).
  * [15] Ee, H.-S. & Agarwal, R.  Tunable metasurface and flat optical zoom lens on a stretchable substrate.  _Nano Letters_ 16, 2818–2823, DOI: 10.1021/acs.nanolett.6b00618 (2016).  PMID: 26986191.
  * [16] Khorasaninejad, M. & Capasso, F.  Metalenses: Versatile multifunctional photonic components.  _Science_ 358, eaam8100, DOI: 10.1126/science.aam8100 (2017).
  * [17] Zhu, B. O. _et al._ Dynamic control of electromagnetic wave propagation with the equivalent principle inspired tunable metasurface.  _Scientific Reports_ 4, 4971, DOI: 10.1038/srep04971 (2014).
  * [18] Cai, T. _et al._ Ultra-thin polarization beam splitter using 2-d transmissive phase gradient metasurface.  _IEEE Transactions on Antennas and Propagation_ 63, 5629–5636, DOI: 10.1109/TAP.2015.2496115 (2015).
  * [19] Yang, J. _et al._ Folded transmitarray antenna with circular polarization based on metasurface.  _IEEE Transactions on Antennas and Propagation_ 69, 806–814, DOI: 10.1109/TAP.2020.3016170 (2021).
  * [20] Yang, Y. _et al._ Dielectric meta-reflectarray for broadband linear polarization conversion and optical vortex generation.  _Nano Letters_ 14, 1394–1399, DOI: 10.1021/nl4044482 (2014).
  * [21] Park, J., Kang, J.-H., Kim, S. J., Liu, X. & Brongersma, M. L.  Dynamic reflection phase and polarization control in metasurfaces.  _Nano Letters_ 17, 407–413, DOI: 10.1021/acs.nanolett.6b04378 (2017).
  * [22] Nayeri, P., Yang, F. & Elsherbeni, A. Z.  Beam-scanning reflectarray antennas: A technical overview and state of the art.  _IEEE Antennas and Propagation Magazine_ 57, 32–47, DOI: 10.1109/MAP.2015.2453883 (2015).
  * [23] Yang, X. _et al._ A broadband high-efficiency reconfigurable reflectarray antenna using mechanically rotational elements.  _IEEE Transactions on Antennas and Propagation_ 65, 3959–3966, DOI: 10.1109/TAP.2017.2708079 (2017).
  * [24] Faenzi, M. _et al._ Metasurface antennas: New models, applications and realizations.  _Scientific Reports_ 9, 10178, DOI: 10.1038/s41598-019-46522-z (2019).
  * [25] Moeini, M. M., Oraizi, H., Amini, A. & Nayyeri, V.  Wide-band beam-scanning by surface wave confinement on leaky wave holograms.  _Scientific Reports_ 9, 13227, DOI: 10.1038/s41598-019-49619-7 (2019).
  * [26] Sievenpiper, D., Colburn, J., Fong, B., Ottusch, J. & Visher, J.  Holographic artificial impedance surfaces for conformal antennas.  In _2005 IEEE Antennas and Propagation Society International Symposium_ , vol. 1B, 256–259 vol. 1B, DOI: 10.1109/APS.2005.1551536 (2005).
  * [27] Fong, B. H., Colburn, J. S., Ottusch, J. J., Visher, J. L. & Sievenpiper, D. F.  Scalar and tensor holographic artificial impedance surfaces.  _IEEE Transactions on Antennas and Propagation_ 58, 3212–3221, DOI: 10.1109/TAP.2010.2055812 (2010).
  * [28] Li, Y. B., Wan, X., Cai, B. G., Cheng, Q. & Cui, T. J.  Frequency-controls of electromagnetic multi-beam scanning by metasurfaces.  _Scientific Reports_ 4, 6921, DOI: 10.1038/srep06921 (2014).
  * [29] Teniou, M., Roussel, H., Capet, N., Piau, G.-P. & Casaletti, M.  Implementation of radiating aperture field distribution using tensorial metasurfaces.  _IEEE Transactions on Antennas and Propagation_ 65, 5895–5907, DOI: 10.1109/TAP.2017.2751652 (2017).
  * [30] Minatti, G. _et al._ Modulated metasurface antennas for space: Synthesis, analysis and realizations.  _IEEE Transactions on Antennas and Propagation_ 63, 1288–1300, DOI: 10.1109/TAP.2014.2377718 (2015).
  * [31] Amini, A., Oraizi, H., Hamedani, M. & Keivaan, A.  Wide-band polarization control of leaky waves on anisotropic holograms.  _Phys. Rev. Applied_ 13, 014038, DOI: 10.1103/PhysRevApplied.13.014038 (2020).
  * [32] Minatti, G., Caminita, F., Martini, E. & Maci, S.  Flat optics for leaky-waves on modulated metasurfaces: Adiabatic floquet-wave analysis.  _IEEE Transactions on Antennas and Propagation_ 64, 3896–3906, DOI: 10.1109/TAP.2016.2590559 (2016).
  * [33] Amini, A. & Oraizi, H.  Adiabatic floquet-wave expansion for the analysis of leaky-wave holograms generating polarized vortex beams.  _Phys. Rev. Applied_ 15, 034049, DOI: 10.1103/PhysRevApplied.15.034049 (2021).
  * [34] González-Ovejero, D. & Maci, S.  Gaussian ring basis functions for the analysis of modulated metasurface antennas.  _IEEE Transactions on Antennas and Propagation_ 63, 3982–3993, DOI: 10.1109/TAP.2015.2442585 (2015).
  * [35] Bodehou, M., González-Ovejero, D., Craeye, C. & Huynen, I.  Method of moments simulation of modulated metasurface antennas with a set of orthogonal entire-domain basis functions.  _IEEE Transactions on Antennas and Propagation_ 67, 1119–1130, DOI: 10.1109/TAP.2018.2880075 (2019).
  * [36] Cassedy, E. S.  Waves guided by a boundary with time-space periodic modulation.  _Proceedings of the Institution of Electrical Engineers_ 112, 269–279, DOI: 10.1049/piee.1965.0040 (1965).
  * [37] Tretyakov, S.  _Analytical Modeling in Applied Electromagnetics_ (Artech House, Boston, 2003).
  * [38] Oliner, A. & Hessel, A.  Guided waves on sinusoidally-modulated reactance surfaces.  _IRE Transactions on Antennas and Propagation_ 7, 201–208, DOI: 10.1109/TAP.1959.1144771 (1959).
  * [39] Oraizi, H., Amini, A., Abdolali, A. & Karimimehr, A. M.  Design of wideband leaky-wave antenna using sinusoidally modulated impedance surface based on the holography theory.  _IEEE Antennas and Wireless Propagation Letters_ 17, 1807–1811, DOI: 10.1109/LAWP.2018.2866712 (2018).
  * [40] Casaletti, M.  Guided waves on scalar and tensorial reactance surfaces modulated by periodic functions: A circuital approach.  _IEEE Access_ 7, 68823–68836, DOI: 10.1109/ACCESS.2019.2919495 (2019).

Figure 1: Surface impedance boundary condition implemented on a grounded
dielectric slab. The propagation direction is assumed to be along the x-axis.
Figure 2: Scalar surface impedance patterns for different time slots. The
circular indicator shows the progression of impedance with time.

(a)

(b)

(c)

Figure 3: Dispersion curves for different pumping frequencies. (a) $f_{p}=$,
(b) $f_{p}=1GHz$, and (c) $f_{p}=2GHz$.

(a)

(b)

Figure 4: (a) Dispersion diagram of proposed hologram in the pure-space
modulation mode ($f_{p}=0\,GHz$). (b) corresponding bounded and radiation
fields at 18 GHz.

(a)

(b)

(c)

(d)

Figure 5: Field distribution for different values of modulation depth. (a)
$M=0.01$, (b) $M=0.1$, (c) $M=0.2$, and (d) $M=0.6$.

(a)

(b)

Figure 6: Dispersion diagram of the proposed spatiotemporally modulated
hologram for the cas of $f_{p}=2\,GHz$ in (a) transmission mode, and (b)
receiving mode.

(a)

(b)

(c)

(d)

Figure 7: Analytical results of the proposed spationtemporally modulated
hologram in transmission mode. (a) total field, (b) fundamental mode, (c)
radiation mode, and (d) (+1)-index mode.

(a)

(b)

(c)

(d)

Figure 8: Analytical results of the proposed spationtemporally modulated
hologram in receiving mode. (a) total field, (b) fundamental mode, (c)
(-1)-index mode, and (d) radiation mode.  Figure 9: Dispersion diagram of
proposed hologram for different values of $f_{p}$.

(a)

(b)

(c)

(d)

Figure 10: Total fields for different values of pumping frequencies. (a)
$f_{p}=0.5\,GHz$, (b) $f_{p}=1\,GHz$, (c) $f_{p}=3\,GHz$, and (d)
$f_{p}=4\,GHz$.


# SMT-based Weighted Model Integration with Structure Awareness

Giuseppe Spallitta University of Trento  Gabriele Masina  University of Trento
Paolo Morettin  KU Leuven  Andrea Passerini  University of Trento  Roberto
Sebastiani  University of Trento

###### Abstract

Weighted Model Integration (WMI) is a popular formalism aimed at unifying
approaches for probabilistic inference in hybrid domains, involving logical
and algebraic constraints. Despite a considerable amount of recent work,
allowing WMI algorithms to scale with the complexity of the hybrid problem is
still a challenge. In this paper we highlight some substantial limitations of
existing state-of-the-art solutions, and develop an algorithm that combines
SMT-based enumeration, an efficient technique in formal verification, with an
effective encoding of the problem structure. This allows our algorithm to
avoid generating redundant models, resulting in substantial computational
savings. An extensive experimental evaluation on both synthetic and real-world
datasets confirms the advantage of the proposed solution over existing
alternatives.

##  1 Introduction

Weighted Model Integration [Belle et al., 2015] recently emerged as a unifying
formalism for probabilistic inference in hybrid domains, characterized by both
continuous and discrete variables and their relationships. The paradigm
extends Weighted Model Counting (WMC) [Chavira and Darwiche, 2008], which is
the task of computing the weighted sum of the set of satisfying assignments of
a propositional formula, to deal with SMT formulas (e.g. [Barrett et al.,
2009]). Whereas WMC can be made extremely efficient by leveraging component
caching techniques [Sang et al., 2004, Bacchus et al., 2009], these strategies
are hard to apply for WMI because of the tight coupling induced by the
arithmetic constraints. Indeed, component caching approaches for WMI are
restricted to fully factorized densities with few dependencies among
continuous variables [Belle et al., 2016]. Another direction specifically
targets acyclic [Zeng and Van den Broeck, 2019, Zeng et al., 2020a] or loopy
[Zeng et al., 2020b] pairwise models.

Exact solutions for more general classes of densities and constraints mainly
leverage advancements in SMT technology or in knowledge compilation (KC)
[Darwiche and Marquis, 2002]. WMI-PA [Morettin et al., 2017, 2019] relies on
SMT-based Predicate Abstraction (PA) [Lahiri et al., 2006] to reduce the
number of models to be generated and integrated over, and was shown to achieve
substantial improvements over previous solutions. However, we show how WMI-PA
has the major drawback of ignoring the structure of the weight function when
pruning away redundant models. This seriously affects its simplification power
when dealing with symmetries in the density. The use of KC for hybrid
probabilistic inference was pioneered by Sanner and Abbasnejad [2012] and
further refined in a series of later works [Kolb et al., 2018, Dos Martires et
al., 2019, Kolb et al., 2020, Feldstein and Belle, 2021]. By compiling a
formula into an algebraic circuit, KC techniques can exploit the structure of
the problem to reduce the size of the resulting circuit, and are at the core
of many state-of-the-art approaches for WMC [Chavira and Darwiche, 2008].
However, even the most recent solutions for WMI [Dos Martires et al., 2019,
Kolb et al., 2020] have serious troubles in dealing with densely coupled
problems, resulting in exponentially large circuits.

In this paper we introduce a novel algorithm for WMI that aims to combine the
best of both worlds, by introducing weight-structure awareness into PA-based
WMI. The main idea is to iteratively build a formula which mimics the
conditional structure of the weight function, so as to drive the SMT-based
enumeration algorithm preventing it from generating redundant models. An
extensive experimental evaluation on synthetic and real-world datasets shows
substantial computational advantages of the proposed solution over existing
alternatives for the most challenging settings.

Our main contributions can be summarized as follows:

  * •

We identify major efficiency issues of existing state-of-the-art WMI
approaches, both PA and knowledge-compilation based.

  * •

We introduce SA-WMI-PA, a novel WMI algorithm that combines PA with weight-
structure awareness.

  * •

We show how SA-WMI-PA achieves substantial computational improvements over
existing solutions in both synthetic and real-world settings.

##  2 Background

###  2.1 SMT and Predicate Abstraction

Satisfiability Modulo Theories (SMT) (see Barrett et al. [2009]) consists in
deciding the satisfiability of first-order formulas over some given theory.
For the context of this paper, we will refer to quantifier-free SMT formulas
over linear real arithmetic ($\mathcal{LRA}$), possibly combined with
uninterpreted function symbols ($\mathcal{LRA}{}\cup\mathcal{EUF}{}$). We
adopt the notation and definitions in Morettin et al. [2019]. We use
$\mathbb{B}\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}\\{\top,\perp\\}$ to
indicate the set of Boolean value, whereas $\mathbb{R}$ indicates the set of
real values. SMT($\mathcal{LRA}$) formulas combines Boolean variables
$A_{i}\in\mathbb{B}$ and $\mathcal{LRA}$ atoms in the form
$(\sum_{i}c_{i}x_{i}\ \bowtie c)$ (where $c_{i}$ are rational values, $x_{i}$
are real variables in $\mathbb{R}$ and $\bowtie$ is one of the standard
algebraic operators $\\{=,\neq,<,>,\leq,\geq\\}$) by using standard Boolean
operators $\\{\neg,\wedge,\vee,\rightarrow,\leftrightarrow\\}$. In
SMT($\mathcal{LRA}\cup\mathcal{EUF}$), $\mathcal{LRA}$ terms can be
interleaved with uninterpreted function symbols. Some shortcuts are provided
to simplify the reading. The formula $(x_{i}\geq l)\wedge(x_{i}\leq u)$ is
shortened into $\llbracket x_{i}\\!\in[l,u]{}\rrbracket$.

Given an SMT formula $\varphi$, a total truth assignment $\mu$ is a function
that maps every atom in $\varphi$ to a truth value in $\mathbb{B}$. A partial
truth assignment maps only a subset of atoms to $\mathbb{B}$.

###  2.2 Weighted Model Integration (WMI)

Let $\mathbf{x}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{x_{1},...,x_{N}\\}\in\mathbb{R}^{N}$ and
$\mathbf{A}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{A_{1},...,A_{M}\\}\in\mathbb{B}^{M}$ for some integers $N$ and
$M$. $\varphi(\mathbf{x},\mathbf{A})$ denotes an SMT($\mathcal{LRA}$) formula
over variables in $\mathbf{A}$ and $\mathbf{x}$ (subgroup of variables are
admissible), while $w{}(\mathbf{x},\mathbf{A})$ denotes a non-negative weight
function s.t. $\mathbb{R}^{N}\times\mathbb{B}^{M}\longmapsto\mathbb{R}^{+}$.
Intuitively, $w{}$ encodes a (possibly unnormalized) density function over
$\mathbf{A}\times\mathbf{x}$. Hereafter $\mu^{\mathbf{A}}$ denotes a truth
assignment on $\mathbf{A}$, $\mu^{\mathcal{LRA}}$ denotes a truth assignment
on the $\mathcal{LRA}$-atoms of $\varphi$,
$\varphi_{[\mu^{\mathbf{A}}]}(\mathbf{x})$ denotes (any formula equivalent to)
the formula obtained from $\varphi$ by substituting every Boolean value
$A_{i}$ with its truth value in $\mu^{\mathbf{A}}$ and propagating the truth
values through Boolean operators, and
$w_{[\mu^{\mathbf{A}}]}(\mathbf{x})\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}w{}(\mathbf{x},\mu^{\mathbf{A}})$ is $w$ restricted to the truth
values of $\mu^{\mathbf{A}}$.

Given a theory
$\mathcal{T}\in\\{{\mathcal{LRA},\mathcal{LRA}\cup\mathcal{EUF}}\\}$, the
nomenclature
$\mathcal{T}\mathcal{T}\mathcal{A}(\varphi)\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{\mu_{1},...,\mu_{j}\\}$ defines the set of
$\mathcal{T}$-consistent total assignments over both propositional and
$\mathcal{T}$ atoms that propositionally satisfy $\varphi$;
$\mathcal{T}\mathcal{A}(\varphi)\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{\mu_{1},...,\mu_{j}\\}$ represents one set of $\mathcal{T}$
partial assignments over both propositional and $\mathcal{T}$ atoms that
propositionally satisfy $\varphi$, s.t. every total assignment in
$\mathcal{T}\mathcal{T}\mathcal{A}(\varphi)$ is a super-assignment of some of
the partial ones in $\mathcal{T}\mathcal{A}(\varphi)$. Given by
$\varphi(\mathbf{x},\mathbf{A})$, with
$\mathcal{T}\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi)$ we mean the set
of all total truth assignment $\mu^{\mathbf{A}}$ on $\mathbf{A}$ s.t.
$\varphi_{[\mu^{\mathbf{A}}]}(\mathbf{x})$ is $\mathcal{T}$-satisfiable, and
by $\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi)$ a set of partial ones
s.t. every total assignment in
$\mathcal{T}\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi)$ is a super-
assignment of some of the partial ones in
$\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi)$.

The Weighted Model Integral of $w$($\mathbf{x}$,$\mathbf{A}$) over
$\varphi(\mathbf{x},\mathbf{A})$ is defined as follows [Morettin et al.,
2019]:

$\displaystyle\ {\sf WMI}(\varphi,w|\mathbf{x},\mathbf{A})\\!\\!\\!$
$\displaystyle\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}$
$\displaystyle\\!\\!\\!\\!\\!\sum_{\mu^{\mathbf{A}}\in\mathbb{B}^{M}}{\sf
WMI_{nb}}(\varphi_{[\mu^{\mathbf{A}}]},w_{[\mu^{\mathbf{A}}]}|\mathbf{x}),$
(1) $\displaystyle=$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\sum_{\mu^{\mathbf{A}}\in\mathcal{T}\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi)}\\!\\!\\!\\!\\!\\!\\!\\!{\sf
WMI_{nb}}(\varphi_{[\mu^{\mathbf{A}}]},w_{[\mu^{\mathbf{A}}]}|\mathbf{x})$ (2)
$\displaystyle{\sf WMI_{nb}}(\varphi,w|\mathbf{x})\ $
$\displaystyle\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}$
$\displaystyle\int_{\varphi(\mathbf{x})}w(\mathbf{x})\ \mathrm{d}\mathbf{x},$
(3) $\displaystyle\\!\\!=\\!\\!$
$\displaystyle\\!\\!\\!\\!\sum_{\mu^{\mathcal{LRA}}\in\mathcal{T}\mathcal{A}(\varphi)}\int_{\mu^{\mathcal{LRA}}}w(\mathbf{x})\
\mathrm{d}\mathbf{x},$ (4)

where the $\mu^{\mathbf{A}}$’s are all total truth assignments on
$\mathbf{A}$, ${\sf WMI_{nb}}(\varphi,w|\mathbf{x})$ is the integral of
$w{}(\mathbf{x})$ over the set $\\{{\mathbf{x}\ |\ \varphi(\mathbf{x})\ is\
true}\\}$ (“nb” means “no-Booleans”).

We call a support of a weight function $w(\mathbf{x},\mathbf{A})$ any subset
of $\mathbb{R}^{N}\times\mathbb{B}^{M}$ out of which
$w(\mathbf{x},\mathbf{A})=0$, and we represent it as a $\mathcal{LRA}$-formula
$\chi(\mathbf{x},\mathbf{A})$. We recall that, consequently,

$\displaystyle{\sf WMI}(\varphi\wedge\chi,w|\mathbf{x},\mathbf{A})={\sf
WMI}(\varphi,w|\mathbf{x},\mathbf{A}).$ (5)

We consider the class of feasibly integrable on $\mathcal{LRA}$ (${\sf
FI}^{\mathcal{LRA}}$) functions $w(\mathbf{x})$, which contain no conditional
component, and for which there exists some procedure able to compute ${\sf
WMI_{nb}}(\mu^{\mathcal{LRA}},w|\mathbf{x})$ for every set of $\mathcal{LRA}$
literals on $\mathbf{x}$. (E.g., polynomials are ${\sf FI}^{\mathcal{LRA}}$.)
Then we call a weight function $w(\mathbf{x},\mathbf{A})$, feasibly integrable
under $\mathcal{LRA}$ conditions (${\sf FIUC}^{\mathcal{LRA}}$) iff it can be
described in terms of a support $\mathcal{LRA}$-formula
$\chi(\mathbf{x},\mathbf{A})$ ($\top$ if not provided), a set
$\mathbf{\Psi}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{{\psi_{i}(\mathbf{x},\mathbf{A})}\\}_{i=1}^{K}$ of
$\mathcal{LRA}$ conditions, in such a way that, for every total truth
assignment $\mu^{\mathbf{\Psi}}$ to $\mathbf{\Psi}$,
$w_{[\mu^{\mathbf{\Psi}}]}$ ($\mathbf{x}$) is total and ${\sf
FI}^{\mathcal{LRA}}$ in the domain given by the values of
$\langle{\mathbf{x},\mathbf{A}}\rangle$ which satisfy
$(\chi\wedge\mu^{\mathbf{\Psi}})_{[\mu^{\mathbf{A}}]}$. ${\sf
FIUC}^{\mathcal{LRA}}$ functions are all the weight functions which can be
described by means of arbitrary combinations of nested if-then-else’s on
conditions in $\mathbf{A}$ and $\mathbf{\Psi}$, s.t. each branch
$\mu^{\mathbf{\Psi}}$ results into a ${\sf FI}^{\mathcal{LRA}}$ weight
function. Each $\mu^{\mathbf{\Psi}}$ describes a portion of the domain of $w$,
inside which $w_{[\mu^{\mathbf{\Psi}}]}$ ($\mathbf{x}$) is ${\sf
FI}^{\mathcal{LRA}}$, and we say that $\mu^{\mathbf{\Psi}}$ identifies
$w_{[\mu^{\mathbf{\Psi}}]}$ in $w$.

In what follows we assume w.l.o.g. that ${\sf FIUC}^{\mathcal{LRA}}$ functions
are described as combinations of constants, variables, standard mathematical
operators $+,-,\cdot,/$ un-conditioned mathematical functions (e.g.,
$exp,sin,...$), conditional expressions in the form $({\sf If}\ \psi_{i}\ {\sf
Then}\ t_{1i}\ {\sf Else}\ {t_{2i}})$ whose conditions $\psi_{i}$ are
$\mathcal{LRA}$ formulas and terms $t_{1i},t_{2i}$ are ${\sf
FIUC}^{\mathcal{LRA}}$.

###  2.3 WMI via Predicate Abstraction

WMI-PA is an efficient WMI algorithm presented in Morettin et al. [2017, 2019]
which exploits SMT-based predicate abstraction. Let $w(\mathbf{x},\mathbf{A})$
be a ${\sf FIUC}^{\mathcal{LRA}}$ function as above. WMI-PA is based on the
fact that

$\displaystyle{\sf WMI}(\varphi,w|\mathbf{x},\mathbf{A})$ $\displaystyle=$
$\displaystyle\sum_{\mu^{\mathbf{A}^{*}}\in\mathcal{T}\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi^{*})}{\sf
WMI_{nb}}(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]},w^{*}_{[\mu^{\mathbf{A}^{*}}]}|\mathbf{x})$
(6) $\displaystyle\varphi^{*}$
$\displaystyle\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}$
$\displaystyle\varphi\wedge\chi\wedge\bigwedge_{k=1}^{K}(B_{k}\leftrightarrow\psi_{k})$
(7)

$\mathbf{A}^{*}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\mathbf{A}\cup\mathbf{B}$ s.t.
$\mathbf{B}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{{B_{1},...,B_{K}}\\}$ are fresh propositional atoms and
$w^{*}(\mathbf{x},\mathbf{A}\cup\mathbf{B})$ is the weight function obtained
by substituting in $w(\mathbf{x},\mathbf{A})$ each condition $\psi_{k}$ with
$B_{k}$.

The pseudocode of WMI-PA is reported in Algorithm 1. First, the problem is
transformed (if needed) by labeling all conditions $\mathbf{\Psi}$ occurring
in $w(\mathbf{x},\mathbf{A})$ with fresh Boolean variables $\mathbf{B}$. After
this preprocessing stage, the set
$\mathcal{M}^{\mathbf{A}^{*}}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\mathcal{T}\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi^{*})$
is computed by invoking SMT-based predicate abstraction Lahiri et al. [2006],
namely $\mathcal{T}\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{*}]}(\mathbf{A}^{*}))$. Then, the algorithm iterates over
each Boolean assignment $\mu^{\mathbf{A}^{*}}$ in
$\mathcal{M}^{\mathbf{A}^{*}}$. $\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]}$ is
simplified by the ${\sf Simplify}$ procedure. Then, if
$\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]}$ is already a conjunction of literals,
the algorithm directly computes its contribution to the volume by calling
${\sf
WMI_{nb}}(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]},w^{*}_{[\mu^{\mathbf{A}^{*}}]}|\mathbf{x})$.
Otherwise, $\mathcal{T}\mathcal{A}(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]})$ is
computed as $\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]}]}(Atoms(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]})))$
to produce partial assignments, and the algorithm iteratively computes
contributions to the volume for each $\mu^{\mathcal{LRA}}$. We refer the
reader to Morettin et al. [2019] for more details.

Notice that in the actual implementation the potentially-large sets
${\mathcal{M}}^{\mathbf{A}^{*}}$ and ${\mathcal{M}}^{\mathcal{LRA}}$ are not
generated explicitly. Rather, their elements are generated, integrated and
then dropped one-by-one, so that to avoid space blowup.

Algorithm 1 WMI-PA($\varphi$, $w$, $\mathbf{x}$, $\mathbf{A}$)

  * 1:

$\langle{\varphi^{*},w^{*},\mathbf{A}^{*}}\rangle\leftarrow{\sf
LabelConditions}(\varphi,w,\mathbf{x},\mathbf{A})$

  * 2:

${\mathcal{M}}^{\mathbf{A}^{*}}\leftarrow\mathcal{T}\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{*}]}(\mathbf{A}^{*}))$

  * 3:

$vol\leftarrow 0$

  * 4:

for $\mu^{\mathbf{A}^{*}}\in{\mathcal{M}}^{\mathbf{A}^{*}}$ do

    * 5:

${\sf Simplify}(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]})$

    * 6:

if ${\sf LiteralConjunction}(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]})$ then

      * 7:

$vol\leftarrow vol+{\sf
WMI_{nb}}(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]},w^{*}_{[\mu^{\mathbf{A}^{*}}]}|\mathbf{x})$

    * 8:

else

      * 9:

${\mathcal{M}}^{\mathcal{LRA}}\leftarrow\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]}]}(Atoms(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]})))$

      * 10:

for $\mu^{\mathcal{LRA}}\in{\mathcal{M}}^{\mathcal{LRA}}$ do

        * 11:

$vol\leftarrow vol+{\sf
WMI_{nb}}(\mu^{\mathcal{LRA}},w^{*}_{[\mu^{\mathbf{A}^{*}}]}|\mathbf{x})$

  * 12:

return $vol$

##  3 Efficiency Issues

###  3.1 Knowledge Compilation

We start our analysis of ${\sf WMI}$ techniques by noticing a major problem
with existing KC approaches for WMI [Dos Martires et al., 2019, Kolb et al.,
2020], in that they tend easily to blow up in space even with simple weight
functions. Consider, e.g., the case in which

$\displaystyle w(\mathbf{x},\mathbf{A})\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\prod_{i=1}^{N}({\sf if}\ \psi_{i}\ {\sf then}\ w_{i1}(\mathbf{x})\
{\sf else}\ w_{i2}(\mathbf{x}))$ (8)

where the $\psi_{i}$s are $\mathcal{LRA}$ conditions on
$\\{{\mathbf{x},\mathbf{A}}\\}$ and the $w_{i1},w_{i2}$ are generic functions
on $\mathbf{x}$. First, the decision diagrams do not interleave arithmetical
and conditional operators, rather they push all the arithmetic operators below
the conditional ones. Thus with (8) the resulting decision diagrams consist of
$2^{N}$ branches on the $\psi_{i}$s, each corresponding to a distinct
unconditioned weight function of the form
$\prod_{i=1}^{N}w_{i{j_{i}}}(\mathbf{x})$ s.t. $j_{i}\in\\{{1,2}\\}$. Second,
the decision diagrams are built on the Boolean abstraction of
$w(\mathbf{x},\mathbf{A})$, s.t. they do not eliminate a priori the useless
branches consisting in $\mathcal{LRA}$-inconsistent combinations of
$\psi_{i}$s, which can be up to exponentially many.

With WMI-PA, instead, the representation of (8) does not grow in size, because
${\sf FIUC}^{\mathcal{LRA}}$ functions allow for interleaving arithmetical and
conditional operators. Also, the SMT-based enumeration algorithm does not
generate $\mathcal{LRA}$-inconsistent assignments on the $\psi_{i}$s. We
stress the fact that (8) is not an artificial scenario: rather, e.g., this is
the case of the real-world logistics problems in Morettin et al. [2019].

###  3.2 WMI-PA

We continue our analysis by noticing a major deficiency also of the WMI-PA
algorithm, that is, it fails to leverage the structure of the weight function
to prune the set of models to integrate over. We illustrate the issue by means
of a simple example (see Figure 1).

Figure 1:  Example highlighting the efficiency issues of the WMI-PA algorithm.
(a) definition of formula $\varphi$ (trivially true) and support $\chi$. (b)
definition of the weight function $w(\mathbf{x},\mathbf{A})$. Round nodes
indicate if-then-else conditions, with true and false cases on the right and
left outgoing edges respectively. Squared nodes indicate ${\sf
FI}^{\mathcal{LRA}}$ weight functions. (c) novel version of the formula
$\varphi^{*}(\mathbf{x},\mathbf{A}\cup\mathbf{B})$ after the application of
the LabelConditions(…) step of WMI-PA. (d) novel version of the weight
function $w^{*}(\mathbf{x},\mathbf{A}\cup\mathbf{B})$, where all
$\mathcal{LRA}$ conditions have been replaced with the fresh Boolean variables
introduced in $\varphi^{*}(\mathbf{x},\mathbf{A}\cup\mathbf{B})$. (e) List of
assignments obtained by WMI-PA on $\mathbf{A}\cup\mathbf{B}$. Notice the
amount of assignments sharing the same ${\sf FI}^{\mathcal{LRA}}$ weight
function.

###### Example 1

Let $\varphi\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}\top$,
$\chi\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}\llbracket
x_{1}\\!\in[0,2]{}\rrbracket\wedge\llbracket x_{2}\\!\in[0,3]{}\rrbracket$
(Figure 1(a)) and let $w(\mathbf{x},\mathbf{A})$ be a tree-structured weight
function defined as in Figure 1(b). To compute ${\sf
WMI}(\varphi\wedge\chi,w|\mathbf{x},\mathbf{A})$, six integrals have to be
computed:

$x_{1}^{2}x_{2}$ on $\llbracket x_{1}\\!\in[1,2]{}\rrbracket\wedge\llbracket
x_{2}\\!\in[1,3]{}\rrbracket$ (if $A_{1}=\top$)
$x_{1}^{3}x_{2}$ on $\llbracket x_{1}\\!\in[1,2]{}\rrbracket\wedge\llbracket
x_{2}\\!\in[0,1]{}\rrbracket$ (if $A_{1}=\top$)
$x_{1}x_{2}^{2}$ on $\llbracket x_{1}\\!\in[0,1]{}\rrbracket\wedge\llbracket
x_{2}\\!\in[2,3]{}\rrbracket$ (if $A_{1}=\top$)
$x_{1}x_{2}^{3}$ on $\llbracket x_{1}\\!\in[0,1]{}\rrbracket\wedge\llbracket
x_{2}\\!\in[0,2]{}\rrbracket$ (if $A_{1}=\top$).
$2x_{1}x_{2}$ on $\llbracket x_{1}\\!\in[0,2]{}\rrbracket\wedge\llbracket
x_{2}\\!\in[0,3]{}\rrbracket$ (if $A_{1}=\bot,A_{2}=\top$)
$3x_{1}x_{2}$ on $\llbracket x_{1}\\!\in[0,2]{}\rrbracket\wedge\llbracket
x_{2}\\!\in[0,3]{}\rrbracket$ (if $A_{1}=\bot,A_{2}=\bot$)

When WMI-PA is used (Algorithm 1), applying LabelConditions(…) we obtain
(Figure 1(c)):

$\displaystyle\varphi^{*}(\mathbf{x},\mathbf{A}\cup\mathbf{B})$
$\displaystyle\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}$
$\displaystyle\varphi\wedge\chi\wedge B_{1}\leftrightarrow(x_{1}\geq 1)\wedge$
$\displaystyle B_{2}\leftrightarrow(x_{2}\geq 1)\wedge
B_{3}\leftrightarrow(x_{2}\geq 2)$

and the weight function $w^{*}(\mathbf{x},\mathbf{A}\cup\mathbf{B})$ shown in
Figure 1(d). Then, by applying $\mathcal{T}\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{*}]}(\mathbf{A}^{*}))$ (row 2) we obtain 24 total
assignments ${\mathcal{M}}^{\mathbf{A}^{*}}{}$ on $\mathbf{A}\cup\mathbf{B}$,
as shown in Figure 1(e).. Notice that WMI-PA uselessly splits into 2 parts the
integrals on $x_{1}^{2}x_{2}$ and $x_{1}x_{2}^{3}$ and into 6 parts the
integral on $2x_{1}x_{2}$ and on $3x_{1}x_{2}$. Also, it repeats the very same
integrals for $\\{{A_{1},A_{2},...}\\}$ and $\\{{A_{1},\neg A_{2},...}\\}$.
$\diamond$

We highlight two facts. First, WMI-PA enumerates total truth assignments on
the Boolean atoms $\mathbf{A}\cup\mathbf{B}$ in
$\mathcal{T}\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi^{*})$ (6) (row 2
in Algorithm 1), assigning also unnecessary values. Second, WMI-PA labels
$\mathcal{LRA}$ conditions in $w$ by means of fresh Boolean atoms $\mathbf{B}$
(row 1 in Algorithm 1). This forces the enumerator to assign all their values
in every assignment, even when not necessary.

The key issue about WMI-PA is that the enumeration of
$\mathcal{T}\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi^{*})$ in (6) and
of $\mathcal{T}\mathcal{A}(\varphi^{*}_{[\mu^{\mathbf{A}^{*}}]})$ in (4) (rows
2 and 11 in Algorithm 1) is not aware of the conditional structure of the
weight function $w$, in particular, it is not aware of the fact that often
partial assignments to the set of conditions in $w^{*}$ (both Boolean and
$\mathcal{LRA}$) are sufficient to identify the value of a ${\sf
FIUC}^{\mathcal{LRA}}$ function (e.g $\\{{A_{1},B_{1},B_{2}}\\}$ suffices to
identify $x_{1}^{2}x_{2}$, or $\\{{\neg A_{1},A_{2}}\\}$ suffices to identify
$2x_{1}x_{2}$), so that it is forced to enumerate all total assignments
extending them (e.g. $\\{{A_{1},A_{2},B_{1},B_{2},B_{3}}\\}$ and
$\\{{A_{1},A_{2},B_{1},B_{2},\neg B_{3}}\\}$).

Thus, to cope with this issue, we need to modify WMI-PA to make it aware of
the conditional structure of $w$.

##  4 Making WMI-PA Weight-Structure Aware

The key idea to prevent total enumeration works as follows. We do not rename
with $\mathbf{B}$ the conditions $\mathbf{\Psi}$ in $w$ and, rather than
enumerating total truth assignments for $\exists\mathbf{x}.\varphi^{*}$ as in
(6)–(7), we enumerate partial assignments for $\exists\mathbf{x}.\varphi^{**}$
where $\varphi^{**}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\varphi\wedge\chi\wedge{\sf sk}(w)$ and ${\sf sk}(w)$ –which we
call the conditional skeleton of $w$– is a $\mathcal{LRA}$ formula s.t.:
$(a)$ its atoms are all and only the conditions in $\mathbf{\Psi}$,
$(b)$ is $\mathcal{LRA}$-valid, so that $\varphi\wedge\chi$ is equivalent to
$\varphi\wedge\chi\wedge{\sf sk}(w)$,
$(c)$ any partial truth value assignment $\mu$ to the conditions
$\mathbf{\Psi}$ which makes ${\sf sk}(w)$ true is such that $w_{[\mu]}$ is
${\sf FI}^{\mathcal{LRA}}$.111E.g., the partial assignment
$\mu\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}\\{{A_{1},(x_{1}\geq
1),(x_{2}\geq 1)}\\}$ in Example 1 is such that $w_{[\mu]}=x_{1}^{2}x_{2}$,
which is ${\sf FI}^{\mathcal{LRA}}$.
Thus, we have that (2) can be rewritten as:

$\displaystyle{\sf WMI}(\varphi,w|\mathbf{x},\mathbf{A})$ $\displaystyle=$
$\displaystyle\sum_{\mu\in\mathcal{T}\mathcal{A}(\exists\mathbf{x}.\varphi^{**})}2^{|\mathbf{A}\setminus\mu|}\cdot{\sf
WMI_{nb}}(\varphi^{**}_{[\mu]},w_{[\mu]}|\mathbf{x})$ (9)
$\displaystyle\varphi^{**}$ $\displaystyle\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}$ $\displaystyle\varphi\wedge\chi\wedge{\sf sk}(w)$ (10)

where $|\mathbf{A}\setminus\mu|$ is the number of Boolean atoms in
$\mathbf{A}$ that are not assigned by $\mu$. Condition $(c)$ guarantees that
${\sf WMI_{nb}}(\varphi^{**}_{[\mu]},w_{[\mu]}|\mathbf{x})$ in (9) can be
directly computed, without further partitioning. The
$2^{|\mathbf{A}\setminus\mu|}$ factor in (9) resembles the fact that, if some
Boolean atom $A_{i}\in\mathbf{A}$ is not assigned in $\mu$, then ${\sf
WMI_{nb}}(\varphi^{**}_{[\mu]},w_{[\mu]}|\mathbf{x})$ should be counted twice
because $\mu$ represents two assignments $\mu\cup\\{{A_{i}}\\}$ and
$\mu\cup\\{{\neg A_{i}}\\}$ which would produce two identical integrals.

Notice that logic-wise ${\sf sk}(w)$ is non-informative because it is a valid
formula. Nevertheless, the role of ${\sf sk}(w)$ is to mimic the structure of
$w$ so that to “make the enumerator aware of the presence of the conditions
$\mathbf{\Psi}{}$”, forcing every assignment $\mu$ to assign truth values also
to these conditions which are necessary to make $w_{[\mu]}$ ${\sf
FI}^{\mathcal{LRA}}$ and hence make ${\sf
WMI_{nb}}(\varphi^{**}_{[\mu]},w_{[\mu]}|\mathbf{x})$ directly computable,
without further partitioning.

An important issue is to avoid ${\sf sk}(w)$ blow up in size. E.g., one could
use as ${\sf sk}(w)$ a formula encoding the conditional structure of an XADDs
or (F)XSDDs, but this may cause a blow up in size, as discussed in Section
3.1.

In order to prevent such problems, we do not generate ${\sf sk}(w)$
explicitly. Rather, we build it as a disjunction of partial assignments over
$\mathbf{\Psi}$ which we enumerate progressively. To this extent, we define
${\sf sk}(w)\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\exists\mathbf{y}.\llbracket y=w\rrbracket{}$ where $\llbracket
y=w\rrbracket$ is a formula on $\mathbf{A}$,$\mathbf{x}$,$\mathbf{y}$ s.t.
$\mathbf{y}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{{y,y_{1},...,y_{k}}\\}$ is a set of fresh variables. Thus,
$\mathcal{T}\mathcal{A}(\exists\mathbf{x}.(\varphi\wedge\chi\wedge\exists\mathbf{y}.\llbracket
y=w\rrbracket{}))$ can be computed as
$\mathcal{T}\mathcal{A}(\exists\mathbf{x}\mathbf{y}.(\varphi\wedge\chi\wedge\llbracket
y=w\rrbracket{}))$ because the $\mathbf{y}$’s do not occur in
$\varphi\wedge\chi$, with no need to generate ${\sf sk}(w)$ explicitly. The
enumeration of
$\mathcal{T}\mathcal{A}(\exists\mathbf{x}\mathbf{y}.(\varphi\wedge\chi\wedge\llbracket
y=w\rrbracket{}))$ is performed by the very same SMT-based procedure used in
Morettin et al. [2019].

$\llbracket y=w\rrbracket$ is obtained by taking $(y=w)$, s.t. $y$ is fresh,
and recursively substituting bottom-up every conditional term $({\sf If}\
\psi_{i}\ {\sf Then}\ t_{i1}\ {\sf Else}\ {t_{i2}})$ in it with a fresh
variable $y_{i}\in\mathbf{y}$, adding the definition of $(y_{i}=({\sf If}\
\psi_{i}\ {\sf Then}\ t_{i1}\ {\sf Else}\ {t_{i2}}))$ as

$\displaystyle(\neg\psi_{i}\vee y_{i}=t_{i1})\wedge(\psi_{i}\vee
y_{i}=t_{i2}).$ (11)

This labeling&rewriting process, which is inspired to labeling CNF-ization
[Tseitin, 1968], guarantees that the size of $\llbracket y=w\rrbracket$ is
linear wrt. that of $w$. E.g., if (8) holds, then $\llbracket y=w\rrbracket$
is $(y=\prod_{i=1}^{N}y_{i})\wedge\bigwedge_{i=1}^{N}((\neg\psi_{i}\vee
y_{i}=w_{i1}(\mathbf{x}))\wedge(\psi_{i}\vee y_{i}=w_{i2}(\mathbf{x})))$.

One problem with the above definition of $\llbracket y=w\rrbracket$ is that it
is not a $\mathcal{LRA}$-formula, because $w$ may include multiplications or
even transcendental functions out of the conditions $\mathbf{\Psi}$ 222The
conditions in $\mathbf{\Psi}$ contain only linear terms by definition., which
makes SMT reasoning over it dramatically hard or even undecidable. We notice,
however, that when computing
$\mathcal{T}\mathcal{A}(\exists\mathbf{x}\mathbf{y}.(\varphi\wedge\chi\wedge\llbracket
y=w\rrbracket{}))$ the arithmetical functions (including operators $+,-,*,/$)
occurring in $w$ out of the conditions $\mathbf{\Psi}$ have no role, since the
only fact that we need to guarantee for the validity of ${\sf sk}(w)$ is that
they are indeed functions, so that $\exists y.(y=f(...))$ is always valid.
333This propagates down on the recursive structure of $\llbracket
y=w\rrbracket$ because, if $y$ does not occur in $\psi$, $\exists y.(y=({\sf
If}\ \psi\ {\sf Then}\ t_{1}\ {\sf Else}\ {t_{2}}))$ is equivalent to $(({\sf
If}\ \psi\ {\sf Then}\ \exists y.(y=t_{1})\ {\sf Else}\ {\exists
y.(y=t_{2})}))$, and $\varphi$ is equivalent to $\exists
y.(\varphi[t|y]\wedge(y=t))$.  (In substance, during the enumeration we are
interested only in the truth values of the conditions $\mathbf{\Psi}$ in $\mu$
which make $w_{[\mu]}$ ${\sf FI}^{\mathcal{LRA}}$, regardless the actual
values of $w_{[\mu]}$). Therefore we can safely substitute condition-less
arithmetical functions (including operators $+,-,\cdot,/$) with some fresh
uninterpreted function symbols, obtaining a
$\mathcal{LRA}\cup\mathcal{EUF}$-formula $\llbracket
y=w\rrbracket_{\mathcal{EUF}}$, which is relatively easy to solve by standard
SMT solvers [Barrett et al., 2009]. It is easy to see that a partial
assignment $\mu$ evaluating $\llbracket y=w\rrbracket$ to true is
$\mathcal{LRA}$-satisfiable iff its corresponding assignment
$\mu_{\mathcal{EUF}}$ is $\mathcal{LRA}\cup\mathcal{EUF}$-satisfiable.444This
boils down to the fact that $y$ occurs only in the top equation and as such it
is free to assume arbitrary values, and that all arithmetic functions are
total in the domain so that, for every possible values of $\mathbf{x}$ a value
for $\mathbf{y}$ always exists iff there exists in the $\mathcal{EUF}$
version. Therefore, we can modify the enumeration procedure into
$\mathcal{T}\mathcal{A}(\exists\mathbf{x}\mathbf{y}.(\varphi\wedge\chi\wedge\llbracket
y=w\rrbracket_{\mathcal{EUF}}{}))$.

Finally, we enforce the fact that the two branches of an if-then-else are
alternative by adding to (11) a mutual-exclusion constraint
$\neg(y_{i}=t_{i1})\vee\neg(y_{i}=t_{i2})$, so that the choice of the function
is univocally associated to the list of decisions on the $\psi_{i}$s. The
procedure producing $\llbracket y=w\rrbracket_{\mathcal{EUF}}$ is described in
detail in Appendix, Algorithm 1.

Figure 2:  Example of bottom-up procedure for computing the relabeling
function $\llbracket y=w\rrbracket_{\mathcal{EUF}}$. (a) Replacement of ${\sf
FI}^{\mathcal{LRA}}$ weight functions (the leaves of the tree, highlighted in
red) with $\mathcal{EUF}$ function symbols (we dropped the dependency on $x$
for compactness). (b-g) Sequence of relabeling steps. At each step, a
conditional term is replaced by a fresh $\mathcal{LRA}$ variable $y_{i}$. The
encoding of the variable in shown in the upper part, while the lower part
shows the weight function with the branch of the conditional term replaced
with $y_{i}$ (highlighted in red). The last step consists in renaming the top
variable as $y$, so that $y=w(\mathbf{x},\mathbf{A}){}$. The relabeling
function $\llbracket y=w\rrbracket_{\mathcal{EUF}}$ is simply the conjunction
of the encodings in the different steps. Figure 3:  Example of structure-aware
enumeration performed by SA-WMI-PA on the problem in Example 1. (a) Generation
of the first assignment. The assignment is on top, while the bottom part shows
the $\llbracket y=w\rrbracket_{\mathcal{EUF}}$ formula. Colors indicate the
progression of the generation, in terms of atoms added (top) and parts of the
formula to be removed as a consequence (bottom). For the sake of simplicity
all atoms until the next atom in $\Psi\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{{A_{1},A_{2},(x_{1}\geq 1),(x_{2}\geq 1),(x_{2}\geq 2)}\\}$ (if
any) are given the same color. (b) Generation of the second assignment. Note
how the $\llbracket y=w\rrbracket_{\mathcal{EUF}}$ formula is enriched with
the blocking clause $\neg A_{1}\vee\neg(x_{1}\geq 1)\vee\neg(x_{2}\geq 1)$
preventing the first assignment to be generated again. (c) Final result of the
enumeration (which contains six assignments in total). The partial assignments
are obtained by restricting the generated assignments on the conditions in
$\Psi$ (and combining them with the atoms of $\chi$, which here are assigned
deterministically.). For each assignment, the corresponding chain of
equivalences of the $\mathbf{y}$s with the identified leaf ${\sf
FI}^{\mathcal{LRA}}$ function is displayed.

###### Example 2

Consider the problem in Example 1. Figure 2 shows the relabeling process
applied to the weight function $w$. The resulting $\llbracket
y=w\rrbracket_{\mathcal{EUF}}$ formula is:

$\begin{array}[]{ll}\phantom{\wedge}\llbracket
y=w\rrbracket(\mathbf{x}\cup\mathbf{y},\mathbf{A}){}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\\ \phantom{\wedge}(\neg A_{1}\vee\neg(x_{1}\geq
1)\vee\neg(x_{2}\geq 1)&\vee\phantom{\neg}(y_{1}=f_{11}(\mathbf{x})))\\\
\wedge(\neg A_{1}\vee\neg(x_{1}\geq 1)\vee\phantom{\neg}(x_{2}\geq
1)&\vee\phantom{\neg}(y_{1}=f_{12}(\mathbf{x})))\\\ \wedge(\neg
A_{1}\vee\neg(x_{1}\geq
1)\vee\neg(y_{1}=f_{11}(\mathbf{x}))&\vee\neg(y_{1}=f_{12}(\mathbf{x})))\\\
\wedge(\neg A_{1}\vee\phantom{\neg}(x_{1}\geq 1)\vee\neg(x_{2}\geq
2)&\vee\phantom{\neg}(y_{2}=f_{21}(\mathbf{x})))\\\ \wedge(\neg
A_{1}\vee\phantom{\neg}(x_{1}\geq 1)\vee\phantom{\neg}(x_{2}\geq
2)&\vee\phantom{\neg}(y_{2}=f_{22}(\mathbf{x})))\\\ \wedge(\neg
A_{1}\vee\phantom{\neg}(x_{1}\geq
1)\vee\neg(y_{2}=f_{21}(\mathbf{x}))&\vee\neg(y_{2}=f_{22}(\mathbf{x})))\\\
\wedge(\neg A_{1}\vee\neg(x_{1}\geq 1)&\vee\phantom{\neg}(y_{3}=y_{1}))\\\
\wedge(\neg A_{1}\vee\phantom{\neg}(x_{1}\geq
1)&\vee\phantom{\neg}(y_{3}=y_{2}))\\\ \wedge(\neg
A_{1}\vee\neg(y_{3}=y_{1})\vee\neg(y_{3}=y_{2}))&\\\
\wedge(\phantom{\neg}A_{1}\vee\neg
A_{2}&\vee\phantom{\neg}(y_{4}=f_{3}(\mathbf{x})))\\\
\wedge(\phantom{\neg}A_{1}\vee\phantom{\neg}A_{2}&\vee\phantom{\neg}(y_{4}=f_{4}(\mathbf{x})))\\\
\wedge(\phantom{\neg}A_{1}\vee\neg(y_{4}=f_{3}(\mathbf{x}))\vee\neg(y_{4}=f_{4}(\mathbf{x})))&\\\
\wedge(\neg A_{1}&\vee\phantom{\neg}(y_{5}=y_{3}))\\\
\wedge(\phantom{\neg}A_{1}&\vee\phantom{\neg}(y_{5}=y_{4}))\\\
\wedge(\neg(y_{5}=y_{3})\vee\neg(y_{5}=y_{4}))&\\\
\wedge(\phantom{\neg}(y=y_{5}))&\\\ \end{array}$

Figure 3 illustrates a possible enumeration process. The algorithm enumerates
partial assignments satisfying $\varphi\wedge\chi\wedge\llbracket
y=w\rrbracket_{\mathcal{EUF}}{}$, restricted on the conditions
$\Psi\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}\\{{A_{1},A_{2},(x_{1}\geq
1),(x_{2}\geq 1),(x_{2}\geq 2)}\\}$, which is equivalent to enumerate
$\mathcal{T}\mathcal{A}(\exists\mathbf{x}\mathbf{y}.(\varphi\wedge\chi\wedge\llbracket
y=w\rrbracket_{\mathcal{EUF}}{}))$. Assuming the enumeration procedure picks
nondeterministic choices following the order of the above set 555Like in
Algorithm 2, we pick Boolean conditions first. and assigning positive values
first, then in the first branch the following satisfying partial assignment is
generated, in order: 666Here nondeterministic choices are underlined. The
atoms in $\chi$ are assigned deterministically.
$\begin{array}[]{l}\chi\cup\\{(y=y_{5}),\underline{A_{1}},(y_{5}=y_{3}),\neg(y_{5}=y_{4}),\underline{(x_{1}\geq
1)},\\\ (y_{3}=y_{1}),\neg(y_{3}=y_{2}),\underline{(x_{2}\geq
1)},(y_{1}=f_{11}(\mathbf{x})),\\\
\neg(y_{1}=f_{12}(\mathbf{x}))\\}.\end{array}$

(Notice that, following the chains of true equalities, we have
$y=y_{5}=y_{3}=y_{1}=f_{11}(\mathbf{x})$.) Then the SMT solver extracts from
it the subset $\\{{A_{1},(x_{1}\geq 1),(x_{2}\geq 1)}\\}$ restricted on the
conditions in $\Psi$. Then the blocking clause $\neg A_{1}\vee\neg(x_{1}\geq
1)\vee\neg(x_{2}\geq 1)$ is added to the formula, which prevents to enumerate
the same subset again. This forces the algorithm to backtrack and generate
$\begin{array}[]{l}\chi\cup\\{(y=y_{5}),\underline{A_{1}},(y_{5}=y_{3}),\neg(y_{5}=y_{4}),(x_{1}\geq
1),\\\ (y_{3}=y_{1}),\neg(y_{3}=y_{2}),{\neg(x_{2}\geq
1)},(y_{1}=f_{12}(\mathbf{x})),\\\
\neg(y_{1}=f_{11}(\mathbf{x}))\\}.\end{array}$ producing the assignment:
$\\{{A_{1},(x_{1}\geq 1),\neg(x_{2}\geq 1)}\\}$. 777We refer the reader to
Lahiri et al. [2006] for more details on the SMT-based enumeration algorithm.

Overall, the algorithm enumerates the following ordered collection of partial
assignments restricted to $Atoms(\varphi\wedge\chi)\cup\mathbf{\Psi}$:
$\begin{array}[]{ll}\chi\cup\\{{\phantom{\neg}A_{1},\phantom{\neg}(x_{1}\geq
1),\phantom{\neg}(x_{2}\geq 1)}\\},&//y=...=f_{11}(\mathbf{x})\\\
\chi\cup\\{{\phantom{\neg}A_{1},\phantom{\neg}(x_{1}\geq 1),\neg(x_{2}\geq
1)}\\},&//y=...=f_{12}(\mathbf{x})\\\
\chi\cup\\{{\phantom{\neg}A_{1},\neg(x_{1}\geq 1),\phantom{\neg}(x_{2}\geq
2)}\\},&//y=...=f_{21}(\mathbf{x})\\\
\chi\cup\\{{\phantom{\neg}A_{1},\neg(x_{1}\geq 1),\neg(x_{2}\geq
2)}\\},&//y=...=f_{22}(\mathbf{x})\\\ \chi\cup\\{{\neg
A_{1},\phantom{\neg}A_{2}}\\},&//y=...=f_{3}(\mathbf{x})\\\ \chi\cup\\{{\neg
A_{1},\neg A_{2}}\\}&//y=...=f_{4}(\mathbf{x})\\\ \end{array}$

which correspond to the six integrals of Example 1. Notice that according to
(2) the first four integrals have to be multiplied by 2, because the partial
assignment $\\{{A_{1}}\\}$ covers two total assignments $\\{{A_{1},A_{2}}\\}$
and $\\{{A_{1},\neg A_{2}}\\}$. Notice also that the disjunction of the six
partial assignments,
$(A_{1}\wedge(x_{1}\geq 1)\wedge(x_{2}\geq 1))\vee...\vee(\neg A_{1}\wedge\neg
A_{2}),$ matches the definition of ${\sf sk}(w)$, which we have computed by
progressive enumeration rather than encoded a priori. $\diamond$

Based on the previous ideas, we develop SA-WMI-PA, a novel “weight structure
aware” variant of WMI-PA. The pseudocode of SA-WMI-PA is reported in Algorithm
2.

As with WMI-PA, we enumerate the assignments in two main steps: in the first
loop (rows 3:-7:) we generate a set ${\mathcal{M}}^{\mathbf{A}^{*}}{}$ of
partial assignments $\mu^{\mathbf{A}^{*}}$ over the Boolean variables
$\mathbf{A}$, s.t. $\varphi^{**}_{[\mu^{\mathbf{A}}]}$ is
$\mathcal{LRA}$-satisfiable and does not contain Boolean variables anymore. In
Example 2 ${\mathcal{M}}^{\mathbf{A}^{*}}{}\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\\{{\\{{A_{1}}\\},\\{{\neg A_{1},A_{2}}\\},\\{{\neg A_{1},\neg
A_{2}}\\}}\\}$. In the second loop (rows 11:-11:), for each
$\mu^{\mathbf{A}^{*}}$ in ${\mathcal{M}}^{\mathbf{A}^{*}}{}$ we enumerate the
set ${\mathcal{M}}^{\mathcal{LRA}}$ of $\mathcal{LRA}$-satisfiable partial
assignments satisfying $\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]}$ (that is, on
$\mathcal{LRA}$ atoms in $Atoms(\varphi\wedge\chi)\cup\mathbf{\Psi}$), we
compute the integral ${\sf
WMI_{nb}}(\mu^{\mathcal{LRA}},w_{[\mu^{\mathbf{A}}*]}|\mathbf{x})$, multiply
it by the $2^{|\mathbf{A}\setminus\mu^{\mathbf{A}^{*}}|}$ factor and add it to
the result. In Example 2, if e.g. $\mu^{\mathbf{A}^{*}}{}=\\{{A_{1}}\\}$,
$\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]}]}(Atoms(\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]})))$
computes the four partial assignments $\\{\chi\cup\\{{(x_{1}\geq 1),(x_{2}\geq
1)}\\},...,\chi\cup\\{{\neg(x_{1}\geq 1),\neg(x_{2}\geq 2)}\\}\\}$.

Algorithm 2 SA-WMI-PA($\varphi$, $w$, $\mathbf{x}$, $\mathbf{A}$)

  * 1:

${\mathcal{M}}^{\mathbf{A}^{*}}\leftarrow\emptyset;vol\leftarrow 0$

  * 2:

$\varphi^{**}\leftarrow\varphi\wedge\chi\wedge\llbracket
y=w\rrbracket_{\mathcal{EUF}}{}$

  * 3:

${\mathcal{M}}^{\mathbf{A}}\leftarrow\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{**}]}(\mathbf{A}))$

  * 4:

for $\mu^{\mathbf{A}}\in{\mathcal{M}}^{\mathbf{A}}$ do

    * 5:

${\sf Simplify}(\varphi^{**}_{[\mu^{\mathbf{A}}]})$

    * 6:

if $\varphi^{**}_{[\mu^{\mathbf{A}}]}$ does not contain Boolean variables then

      * 7:

${\mathcal{M}}^{\mathbf{A}^{*}}\leftarrow{\mathcal{M}}^{\mathbf{A}^{*}}\cup\\{{\mu^{\mathbf{A}}}\\}$

    * 8:

else

      * 9:

for $\mu^{\mathbf{A}}_{residual}\in\mathcal{T}\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{**}_{[\mu^{\mathbf{A}}]}]}(\mathbf{A}))$ do

        * 10:

${\mathcal{M}}^{\mathbf{A}^{*}}\leftarrow{\mathcal{M}}^{\mathbf{A}^{*}}\cup\\{{\mu^{\mathbf{A}}\wedge\mu^{\mathbf{A}}_{residual}}\\}$

  * 11:

for $\mu^{\mathbf{A}^{*}}\in{\mathcal{M}}^{\mathbf{A}^{*}}$ do

    * 12:

$k\leftarrow|\mathbf{A}\setminus\mu^{\mathbf{A}^{*}}|$

    * 13:

${\sf Simplify}(\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]})$

    * 14:

if ${\sf LiteralConjunction}(\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]})$ then

      * 15:

$vol\leftarrow vol+2^{k}\cdot{\sf
WMI_{nb}}(\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]},w_{[\mu^{\mathbf{A}}*]}|\mathbf{x})$

    * 16:

else

      * 17:

${\mathcal{M}}^{\mathcal{LRA}}\leftarrow\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]}]}(Atoms(\varphi^{**}_{[\mu^{\mathbf{A}^{*}}]})))$

      * 18:

for $\mu^{\mathcal{LRA}}\in{\mathcal{M}}^{\mathcal{LRA}}$ do

        * 19:

$vol\leftarrow vol+2^{k}\cdot{\sf
WMI_{nb}}(\mu^{\mathcal{LRA}},w_{[\mu^{\mathbf{A}}*]}|\mathbf{x})$

  * 20:

return $vol$

In detail, in row 2: we extend $\varphi\wedge\chi$ with $\llbracket
y=w\rrbracket_{\mathcal{EUF}}{}$ to provide structure awareness. (We recall
that, unlike with WMI-PA, we do not label $\mathcal{LRA}$ conditions with
fresh Boolean variables $\mathbf{B}{}$.) Next, in row 3: we perform
$\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{**}_{[\mu^{\mathbf{A}}]}]}(\mathbf{A}))$ to obtain a set
${\mathcal{M}}^{\mathbf{A}}{}$ of partial assignments restricted on Boolean
atoms $\mathbf{A}$. Then, for each assignment
$\mu^{\mathbf{A}}\in{\mathcal{M}}^{\mathbf{A}}$ we build the (simplified)
residual $\varphi^{**}_{[\mu^{\mathbf{A}}]}$. Since $\mu^{\mathbf{A}}$ is
partial, $\varphi^{**}_{[\mu^{\mathbf{A}}]}$ is not guaranteed to be free of
Boolean variables $\mathbf{A}$, as shown in Example 3. If this is the case, we
simply add $\mu^{\mathbf{A}}$ to ${\mathcal{M}}^{\mathbf{A}^{*}}$, otherwise
we invoke $\mathcal{T}\mathcal{T}\mathcal{A}({\sf
PredAbs}_{[\varphi^{**}_{[\mu^{\mathbf{A}}]}]}(\mathbf{A}))$ to assign the
remaining variables and conjoin each assignment $\mu^{\mathbf{A}}_{residual}$
to $\mu^{\mathbf{A}}$, ensuring that the residual now contains only
$\mathcal{LRA}$ atoms (rows 4:\- 4:). The second loop (rows 11:-11:) resembles
the main loop in WMI-PA, with the only relevant difference that, since
$\mu^{\mathbf{A}^{*}}$ is partial, the integral is multiplied by a
$2^{|\mathbf{A}\setminus\mu^{\mathbf{A}^{*}}|}$ factor.

Notice that in general the assignments $\mu^{\mathbf{A}^{*}}$ are partial even
if the steps in rows 9:-10: are executed; the set of residual Boolean
variables in $\varphi^{**}_{[\mu^{\mathbf{A}}]}$ are a (possibly much smaller)
subset of $\mathbf{A}\setminus\mu^{\mathbf{A}}$ because some of them do not
occur anymore in $\varphi^{**}_{[\mu^{\mathbf{A}}]}$ after the simplification,
as shown in the following example.

###### Example 3

Let $\varphi\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}(A_{1}\vee A_{2}\vee
A_{3})\wedge(\neg A_{1}\vee A_{2}\vee(x\geq 1))\wedge(\neg A_{2}\vee(x\geq
2))\wedge(\neg A_{3}\vee(x\leq 3))$, $\chi\stackrel{{\scriptstyle\text{\tiny
def}}}{{=}}\llbracket x_{1}\\!\in[0,4]{}\rrbracket$ and
$w(\mathbf{x},\mathbf{A})\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}1.0$.
Suppose $\mathcal{T}\mathcal{A}({\sf PredAbs}_{[\varphi^{**}]}(\mathbf{A}))$
finds the partial assignment $\\{{(x\geq 0),(x\leq 4),A_{2},(x\geq 1),(x\geq
2),(x\leq 3)}\\}$, whose projected version is
$\mu^{\mathbf{A}}\stackrel{{\scriptstyle\text{\tiny def}}}{{=}}\\{{A_{2}}\\}$
(row 3:). Then $\varphi^{**}_{[\mu^{\mathbf{A}}]}$ reduces to $(\neg
A_{3}\vee(x\leq 3))$, so that ${\mathcal{M}}^{\mathbf{A}^{*}}$ is
$\\{{\\{{A_{2},A_{3}}\\},\\{{A_{2},\neg A_{3}}\\}}\\}$, avoiding branching on
$A_{1}$. $\diamond$

We stress the fact that in our actual implementation, like with that of WMI-
PA, the potentially-large sets ${\mathcal{M}}^{\mathbf{A}^{*}}$ and
${\mathcal{M}}^{\mathcal{LRA}}$ are not generated explicitly. Rather, their
elements are generated, integrated and then dropped one-by-one, so that to
avoid space blowup.

We highlight two main differences wrt. WMI-PA. First, unlike with WMI-PA, the
generated assignments $\mu^{\mathbf{A}}$ on $\mathbf{A}$ are partial, each
representing $2^{|\mathbf{A}\setminus\mu^{\mathbf{A}}|}$ total ones. Second,
the assignments on (non-Boolean) conditions $\mathbf{\Psi}$ inside the
$\mu^{\mathcal{LRA}}$s are also partial, whereas with WMI-PA the assignments
to the $\mathbf{B}$s are total. This may drastically reduce the number of
integrals to compute, as empirically demonstrated in the next section.

##  5 Experimental Evaluation

The novel algorithm SA-WMI-PA is compared to the original WMI-PA algorithm
[Morettin et al., 2019], and the WMI solvers based on KC: XADD [Kolb et al.,
2018], XSDD and FXSDD [Kolb et al., 2020]. Each of these methods is called
from the Python framework pywmi [Kolb et al., 2019]. For both WMI-PA and SA-
WMI-PA, we use MathSAT for SMT enumeration and LattE Integrale for computing
integrals. For the KC algorithms we use PSiPSI [Gehr et al., 2016] as symbolic
computer algebra backend. All experiments are performed on an Intel Xeon Gold
6238R @ 2.20GHz 28 Core machine with 128GB of ram and running Ubuntu Linux
20.04. The code of SA-WMI-PA is freely available at https://github.com/unitn-
sml/wmi-pa.
For improved readability, in both the experiments we report runtime using
cactus plots, i.e. the single problem instances are increasingly sorted by
runtime for each algorithm separately. We highlight how, by construction,
problem instances of the same tick of the x-axis are not guaranteed to be the
same for different algorithms. Steeper slopes of an algorithm curve means less
efficiency.

|
---|---
Figure 4: Cactus plots reporting execution time for all methods on the
synthetic experiments (left); number of integrals (right) for WMI-PA and SA-
WMI-PA.

###  5.1 Synthetic experiments

We first evaluate our algorithm on random formulas and weights, following the
experimental protocol of Morettin et al. [2019]. We define two recursive
procedures to generate $\mathcal{LRA}$ formulae and weight functions with
respect to a positive integer number $D$, named depth:

$\displaystyle{\sf r}_{\varphi}(D)=$
$\displaystyle\begin{cases}\bigoplus_{q=1}^{Q}{\sf r}_{\varphi}(D-1)\hskip
65.44142pt\text{ if }D>0\\\ [\neg]A^{\mathbb{B}/\mathbb{R}}\hskip
101.00737pt\text{ otherwise}\\\ \end{cases}$ $\displaystyle{\sf r}_{w}(D)=$
$\displaystyle\begin{cases}\begin{cases}\llbracket{\sf If}\ {\sf
r}_{\varphi}(D)\ {\sf Then}\ {\sf r}_{w}(D-1)\ {\sf Else}\ {\sf
r}_{w}(D-1){}\rrbracket&\\\ \text{or}\\\ {\sf r}_{w}(D-1)\otimes{\sf
r}_{w}(D-1)\ \ \ \ \ \ \ \ \ \ \ \text{ if }D>0\\\ \end{cases}\\\ P_{\sf
random}(\mathbf{x})\hskip 91.04872pt\text{ otherwise}\end{cases}$

where $\bigoplus\in\\{\bigvee,\bigwedge,\neg\bigvee,\neg\bigwedge\\}$,
$\otimes\in\\{+,\cdot\\}$ and $P_{\sf random}(\mathbf{x})$ is a random
polynomial function. Using these procedure we generate instances of synthetic
problems:

$\displaystyle\chi(\mathbf{x},\mathbf{A})$ $\displaystyle=$ $\displaystyle{\sf
r}_{\varphi}(D)\wedge\bigwedge_{x\in\mathbf{x}}\llbracket
x\\!\in[l_{x},u_{x}]{}\rrbracket$ $\displaystyle w(\mathbf{x},\mathbf{A})$
$\displaystyle=$ $\displaystyle{\sf r}_{w}(D)$ $\displaystyle\varphi_{{\sf
query}}(\mathbf{x},\mathbf{A})$ $\displaystyle=$ $\displaystyle{\sf
r}_{\varphi}(D)$

where $l_{x},u_{x}$ are real numbers such that $\forall x.(l_{x}<u_{x})$.

In contrast with the benchmarks used in other recent works [Kolb et al., 2018,
2020], the procedure is not strongly biased towards the generation of problems
with structural regularities, offering a more neutral perspective on how the
different techniques are expected to perform in the wild. The generated
synthetic benchmark contains problems where the number of both Boolean and
real variables is set to 3, while the depth of weights functions fits in the
range $[4,7]$. Timeout was set to 3600 seconds, similarly to what has been
done in previous works.

In this settings, the approaches based on a SMT oracle clearly outperform
those based on KC (Fig. 4 left). In addition, SA-WMI-PA greatly improves over
WMI-PA, thanks to a drastic reduction in the number of integrals computed
(Fig. 4 right), with the advantage of our approach getting more evident when
the weight functions are deeper.

|
---|---
|
Figure 5: Cactus plots representing average query execution times and standard
deviation in seconds on a set of DET problems with
$H\in\\{0.25,0.5,0.75,1\\}$.

###  5.2 Density Estimation Trees

We explore the use of WMI solvers for marginal inference in real-world
probabilistic models. In particular, we considered _Density Estimation Trees_
(DETs) [Ram and Gray, 2011], hybrid density estimators encoding piecewise
constant distributions. Having only univariate conditions in the internal
nodes, DETs natively support tractable inference when single variables are
considered. Answering queries like $Pr(X\leq Y)$ requires instead
marginalizing over an oblique constraint, which is reduced to WMI: $Pr(X\leq
Y)=\frac{WMI((X\leq Y)\wedge\chi_{DET},w_{DET})}{WMI(\chi_{DET},w_{DET})}$.
Being able to address this type of queries is crucial to apply WMI-based
inference to e.g. probabilistic formal verification tasks, involving
constraints that the system should satisfy with high probability.

We considered a selection of hybrid datasets from the UCI repository [Dua and
Graff, 2017], reported in Table 1 in the Appendix. Following the approach of
Morettin et al. [2020], discrete numerical features were relaxed into
continuous variables, while $n$-ary categorical features are one-hot encoded
with $n$ binary variables.

After learning a DET on each dataset, we generated a benchmark of increasingly
complex queries, 5 for each dataset, involving a ratio $H\in[0,1]$ of the
continuous variables. More specifically, the queries are linear inequalities
involving a number of variables $k=max(1,\lfloor H\cdot|\mathbf{x}|\rfloor)$.
Figure 5 depicts the runtime of the algorithms for
$H\in\\{0.25,0.5,0.75,1\\}$. Timeout was set to 1200 seconds. KC approaches
have an edge for the simplest cases ($H\leq 0.5$) in which substantial
factorization of the integrals is possible. Contrarily to many other
probabilistic models, which are akin to the case in Section 3.1, DETs are
well-suited for KC-based inference, due to the absence of arithmetic
operations in the internal nodes. When the coupling between variables
increases, however, the advantage of decomposing the integrals is overweight
by the combinatorial reasoning capabilities of SA-WMI-PA. We remark that SA-
WMI-PA is agnostic of the underlying integration procedure, and thus in
principle it could also incorporate a symbolic integration component.

##  6 Conclusion

We presented the first SMT-based algorithm for WMI that is aware of the
structure of the weight function. This is particularly beneficial when the
piecewise density defined on top of SMT($\mathcal{LRA}$) constraints is deep,
as it is often the case with densities learned from data. We evaluated our
algorithmic ideas on both synthetic and real-world problems, obtaining state-
of-the-art results in many settings of practical interest. Providing
unprecedented scalability in the evaluation of complex probabilistic SMT
problems, this contribution directly impacts the use of WMI for the
probabilistic verification of systems. While the improvements described in
this work drastically reduce the number of integrations required to compute a
weighted model integral, the integration itself remains an obstacle to
scalability. Combining the fast combinatorial reasoning offered by SMT solvers
with symbolic integration is a promising research direction.

###### Acknowledgements.

This research was partially supported by TAILOR, a project funded by EU
Horizon 2020 research and innovation programme under GA No 952215 and from the
European Research Council (ERC) under the European Union’s Horizon 2020
research and innovation programme (grant agreement No. [694980] SYNTH:
Synthesising Inductive Data Models). The authors would like to thank Pedro
Zuidberg Dos Martires for his help with XSDDs and FXSDDs.

## References

  * Bacchus et al. [2009] Fahiem Bacchus, Shannon Dalmao, and Toniann Pitassi.  Solving #SAT and Bayesian inference with backtracking search.  _Journal of Artificial Intelligence Research_ , 34(1):391–442, 2009.
  * Barrett et al. [2009] C. W. Barrett, R. Sebastiani, S. A. Seshia, and C. Tinelli.  Satisfiability Modulo Theories.  In _Handbook of Satisfiability_ , chapter 26, pages 825–885. IOS Press, 2009.
  * Belle et al. [2015] V. Belle, A. Passerini, and G. Van den Broeck.  Probabilistic inference in hybrid domains by weighted model integration.  In _Proceedings of 24th International Joint Conference on Artificial Intelligence (IJCAI)_ , 2015.
  * Belle et al. [2016] Vaishak Belle, Guy Van den Broeck, and Andrea Passerini.  Component caching in hybrid domains with piecewise polynomial densities.  In _AAAI_ , 2016.
  * Chavira and Darwiche [2008] Mark Chavira and Adnan Darwiche.  On probabilistic inference by weighted model counting.  _Artificial Intelligence_ , 172(6-7):772–799, 2008.  ISSN 0004-3702.  10.1016/j.artint.2007.11.002.
  * Darwiche and Marquis [2002] A. Darwiche and P. Marquis.  A knowledge compilation map.  _JAIR_ , 2002.
  * Dos Martires et al. [2019] Pedro Zuidberg Dos Martires, Anton Dries, and Luc De Raedt.  Exact and approximate weighted model integration with probability density functions using knowledge compilation.  In _Proceedings of the AAAI Conference on Artificial Intelligence_ , volume 33, pages 7825–7833, 2019.
  * Dua and Graff [2017] Dheeru Dua and Casey Graff.  UCI machine learning repository, 2017.  URL http://archive.ics.uci.edu/ml.
  * Feldstein and Belle [2021] Jonathan Feldstein and Vaishak Belle.  Lifted reasoning meets weighted model integration.  In _Uncertainty in Artificial Intelligence_ , pages 322–332. PMLR, 2021.
  * Gehr et al. [2016] Timon Gehr, Sasa Misailovic, and Martin Vechev.  Psi: Exact symbolic inference for probabilistic programs.  In _International Conference on Computer Aided Verification_ , pages 62–83. Springer, 2016.
  * Kolb et al. [2018] S. Kolb, M. Mladenov, S. Sanner, V. Belle, and K. Kersting.  Efficient symbolic integration for probabilistic inference.  In _IJCAI_ , 2018.
  * Kolb et al. [2019] Samuel Kolb, Paolo Morettin, Pedro Zuidberg Dos Martires, Francesco Sommavilla, Andrea Passerini, Roberto Sebastiani, and Luc De Raedt.  The pywmi framework and toolbox for probabilistic inference using weighted model integration.  In _Proceedings of the Twenty-Eighth International Joint Conference on Artificial Intelligence, IJCAI-19_ , pages 6530–6532. International Joint Conferences on Artificial Intelligence Organization, 7 2019\.
  * Kolb et al. [2020] Samuel Kolb, Pedro Zuidberg Dos Martires, and Luc De Raedt.  How to exploit structure while solving weighted model integration problems.  In _Uncertainty in Artificial Intelligence_ , pages 744–754. PMLR, 2020.
  * Lahiri et al. [2006] S. K. Lahiri, R. Nieuwenhuis, and A. Oliveras.  SMT techniques for fast predicate abstraction.  In _CAV_ , 2006.
  * Morettin et al. [2017] Paolo Morettin, Andrea Passerini, and Roberto Sebastiani.  Efficient weighted model integration via SMT-based predicate abstraction.  In _Proceedings of the Twenty-Sixth International Joint Conference on Artificial Intelligence, IJCAI-17_ , pages 720–728, 2017.
  * Morettin et al. [2019] Paolo Morettin, Andrea Passerini, and Roberto Sebastiani.  Advanced smt techniques for weighted model integration.  _Artificial Intelligence_ , 275, 04 2019.
  * Morettin et al. [2020] Paolo Morettin, Samuel Kolb, Stefano Teso, and Andrea Passerini.  Learning weighted model integration distributions.  In _Proceedings of the AAAI Conference on Artificial Intelligence_ , volume 34, pages 5224–5231, 2020.
  * Ram and Gray [2011] Parikshit Ram and Alexander G Gray.  Density estimation trees.  In _Proceedings of the 17th ACM SIGKDD international conference on Knowledge discovery and data mining_ , pages 627–635. ACM, 2011.
  * Sang et al. [2004] Tian Sang, Fahiem Bacchus, Paul Beame, Henry A. Kautz, and Toniann Pitassi.  Combining component caching and clause learning for effective model counting.  In _SAT_ , 2004.
  * Sanner and Abbasnejad [2012] S. Sanner and E. Abbasnejad.  Symbolic variable elimination for discrete and continuous graphical models.  In _AAAI_ , 2012.
  * Tseitin [1968] G.S. Tseitin.  _On the complexity of derivations in the propositional calculus_ , page 115–125.  Consultants Bureau, 1968.  Part II.
  * Zeng and Van den Broeck [2019] Z. Zeng and G. Van den Broeck.  Efficient search-based weighted model integration.  In _UAI_ , 2019.
  * Zeng et al. [2020a] Z. Zeng, P. Morettin, F. Yan, A. Vergari, and G. Van den Broeck.  Scaling up hybrid probabilistic inference with logical and arithmetic constraints via message passing.  In _ICML_ , 2020a.
  * Zeng et al. [2020b] Z. Zeng, P. Morettin, F. Yan, A. Vergari, and G. Van den Broeck.  Probabilistic inference with algebraic constraints: Theoretical limits and practical approximations.  In _NeurIPS_ , 2020b.