Consider a two phase direct-search scheme applied to problem (<ref>) using the infeasibility measure (<ref>), and suppose that this
measure is Lipschitz continuous.
Let $\K$ be a sequence of unsuccessful iterations, and let $\{\xb_k\}_{k \in \K}$ be a refining subsequence with refined point $\xb^*$.
For any refining direction $\db$, the following results hold:

* If $h(\xb^*) > 0$, then $h^\circ(\xb^*, \db) \geq 0$.

* If $h(\xb_k) = 0$ for $k$ large enough and $\db \in T^H_{\Omega}(\xb^*)$ then $\xb^*$ is feasible and $f(\xb^*; \db) \geq 0$.

* Assume
$h^\circ(\tilde{\xb}^*; \tilde{\db}) \leq -\varepsilon$  for  $\tilde{\xb}^*$ in a neighborhood of $\xb^*$ intersected with $\mathbb{R}^n \setminus \Omega $ and $\tilde{\db}$ in a neighborhood of $\db$. Then if $h(\xb_k) \rightarrow 0$ we have $h(\xb_k) = 0$ for $k$ large enough, i.e., a feasible point is found by the algorithm.

The first two points directly follow from the results of Theorem <ref> and Theorem <ref>, respectively.
We prove the third point by contradiction. Suppose that $h(\xb_k) > 0$ for every $k$. By the generalized mean value theorem, we get

\begin{equation} \label{eq:meant}
h(\xb_k + \alpha_k \db_k) - h(\xb_k) \leq \alpha_k \sup_{t \in [0, 1]} h^\circ(\xb_k + t \alpha_k; \db_k) \leq -\varepsilon \alpha_k = \Theta(\alpha_k)  \, .
\end{equation}

But since $\rho(\alpha_k) = o(\alpha_k)$ (<ref>) implies that $h(\xb_k + \alpha_k \db_k) - h(\xb_k) < -\rho(\alpha_k)$ for $k \in K$ large enough. This contradicts the assumption that the steps in $K$ are unsuccessful.

Note that a similar result can be obtained using (<ref>) in lieu of (<ref>), assuming differentiability of the function $\bm{g}$ and a weak linear independence condition <cit.>. In both cases, the results illustrate a shortcoming of the approach already mentioned for extreme barrier techniques, namely that the information
at infeasible points is not used. In practice, however, improving the objective function, even while significantly worsening infeasibility, can prove beneficial. In the next three sections, we describe
algorithms based on this paradigm.

§.§.§ Filter-based methods

Filter based approaches roughly deals with the constraints using a bilevel optimization approach, with a set of so called "incumbent" points $S_k$ updated at each iteration, and stored together with the corresponding objective values and infeasibility measures $\{(f(\xb_k^{(i)}), h(\xb_k^{(i)}))\}_{i \in S_k}$. A point $\xb$ is said to dominate a point $\yb$, written $\xb \succ \yb$, if $f(\xb) \geq f(\yb)$ and $h(\xb) \geq h(\yb)$, with one of the two inequalities strict. At every iteration $S_k$ is updated adding some suitably chosen points and discarding the points dominated by some of the newly added points. It is convenient to write $S_k$ as the disjoint union of the feasible incumbent points $I_k$ and of the unfeasible incumbent points $U_k$. The first instance of filter based direct-search method was presented by Audet and Dennis [Audet and Dennis, Jr., 2004], but it was not proved to converge to stationary points for the original problem. Dennis et al. [Dennis, Jr. et al., 2004] described another filter based method was proved to converge to stationary points  assuming however differentiability of the constraints with approximate gradients available, as well as of the objective at limit points. We focus on a third approach
due to Audet and Dennis [Audet and Dennis, Jr., 2009, Audet et al., 2010] called progressive barrier, that combine MADS with a filter technique. A simplified version of this method based on using Algorithm <ref> in lieu of MADS is described in Algorithm <ref>. We choose Algorithm <ref> with sufficient decrease fo simplicity, but note that using Algorithm <ref>
would lead to similar convergence results. The method uses relations of the form $\xb \prec_{\varepsilon} \yb$, which correspond $\xb \prec \yb$ and $f(\xb) < f(\yb) - \varepsilon$ or $h(\xb) < h(\yb)- \varepsilon$, i.e., one of the two inequalities defining majorization is satisfied within a margin of $\varepsilon$.

Progressive barrier direct-search scheme

[1]

Inputs:  			A starting set of incumbents $S_0$, an initial stepsize $\alpha_0> 0$, a parameter
for stepsize contraction $0 < \theta < 1$, a parameter for stepsize expansion
$\gamma \geq 1$, forcing function $\rho: \Rplus \rightarrow \Rplus$, an initial unfeasibility threshold $h_0 > 0$.
Select a set $\mathbb{D}_k$ of poll directions, an unfeasible incumbent poll center $\xb_k^U$ and, if $I_k \neq \emptyset$, a feasible incumbent poll center $\xb_k^I$.

$\yb + \sba \prec_{\rho(\alpha_k)}  \yb$ holds for $\sba \in \{\alpha \db \ | \ \db \in \bbD_k \}$, and $\yb \in \{\xb_k^U, \xb_k^I\}$
Declare the iteration as successful, add $\yb + \sba$ to $S_k$, set $h_{k + 1} = h_k$, and $\alpha_{k + 1} = \gamma \alpha_k$.
Go to step 12.
$h(\xk^U + \bm{s}) < h(\xk^U) - \rho(\alpha_k)$ and  $f(\xk^U + \bm{s}) \leq f(\xk^I)$
hold for some $\bm{s}\in \accolade{\alpha_k\bm{d}:\bm{d} \in \mathbb{D}_k}$
Declare the iteration as improving, add $\xb_k^U + \bm{s}$ to $S_k$, and lower the acceptance threshold by
choosing $h_{k + 1} \in [\bar{h}_k, h_{k} - \rho(\alpha_k)]$. Set $\alpha_{k + 1} = \alpha_k$.  Go to step 12.

Declare the iteration unsuccessful, set $S_{k + 1} = S_k$, $h_{k + 1} = h_k$ and $\alpha_{k + 1} = \theta \alpha_k$.
Delete all the dominated points from $S_{k + 1}$.

The key result for proving convergence of Algorithm <ref> is given below.

Consider a sequence generated by Algorithm <ref>. Assume that, for some $\bar{f} > 0$, $f(\yb) \leq \bar{f}$ for every $\yb \in S_k$. Then we have $\lim \alpha_k = 0$.

Assume, by contradiction, that we have an infinite number of improving or successful steps with $\alpha_k ~>~\varepsilon$, for some $\varepsilon > 0$. First, since $\{h_k\}$ is non increasing and strictly decreases by $\rho(\alpha_k)$ at every improving step, there cannot be an infinite number of improving steps with $\alpha_k > \varepsilon$. Thus we must necessarily have an infinite number of successful steps with $\alpha_k \geq \varepsilon$. Now, on the one hand $ \bar{S} = \cup_k \{ (\lfloor f(\yb)/\rho(\varepsilon) \rfloor, \lfloor h(\yb)/ \rho(\varepsilon) \rfloor) \}_{\yb \in S_k}$ must be bounded because it is a subset of $[0:\lfloor \bar{f}/\rho(\varepsilon) \rfloor] \times [0:\lfloor h_0/\rho(\varepsilon) \rfloor]$. On the other hand, it is easy to check that due to the $\rho(\alpha_k)$ majorization condition at every successful iteration with $\alpha_k \geq \varepsilon$ a new point is added to $\bar{S}$, a contradiction.

Thanks to Proposition <ref>, the results of Theorem <ref> can be applied directly to Algorithm <ref>.

§.§.§ Restoration based methods

Restoration based approaches rely on alternating steps to improve the objective with a restoration
phase focusing on the improvement of the unfeasibility measure $h(\xb)$. Bueno et al.[Bueno et al., 2013] proposed such an approach to deal with differentiable constraints, where a gradient-based restoration phase is alternated with a direct-search phase applied on a problem constrained with linear approximations of the constraints. Our result below is based on the method by Gratton and
Vicente [Gratton and Vicente, 2014], where the following two conditions to determine success of the algorithm at iteration $k$:

\begin{equation}
\label{eq:restoration}
h(\xb_k) - h(\xb_k + \alpha_k \db_k) > \rho(\alpha_k) \quad \text{and} \quad h(\xb_k) > C \rho(\alpha_k)
\end{equation}

for some $C > 1$, and

\begin{equation}
\label{eq:improvement}
M(\xb_k, \mu) - M(\xb_k + \alpha_k \db_k, \mu) > \rho(\alpha_k)
\end{equation}

for some $\mu > 0$, for $M(\xb, \mu) = f(\xb) + \mu h(\xb)$ referred to as merit function. Within the restoration phase, only (<ref>) is considered. While the methodological and theoretical results presented by Gratton and Vicente [Gratton and Vicente, 2014] are beyond the scope of this survey, we show in the next proposition how a convergence result analogous to the one of two phase methods can be obtained by considering suitable sequences of unsuccessful steps.

Assume $\{\xb_k\}_{k \in K}$ is a converging subsequence of unsuccessful iterates where neither (<ref>) nor (<ref>) are
satisfied in tentative points, and let $\db$ be a refined direction for this subsequence. Then for $\xb^*$ limit point of $\{\xb_k\}$ we have
$h^\circ(\xb^*, \db) \geq 0$ and $f^\circ(\xb^*, \db) \geq 0$ whenever $h^\circ(\xb^*, \db) = 0$. In particular,
if $\cup_{k \in K} D_k$ is dense then $\xb^*$ is Clarke stationary for $h$. If furthermore $h$ is a suitably defined infeasibility
measure and $h(\xb^*) = 0$, then $\xb^*$ is a KKT stationary point for the problem.

We get $h^\circ(\xb^*; \db) \geq 0$ and $M^\circ(\xb^*, \mu; \db) \geq 0$  with the usual arguments from the assumption that (<ref>) and (<ref>) respectively are not satisfied for $k \in K$. If $h^\circ(\xb^*; \db) = 0$ then

\begin{equation}
f^\circ(\xb^*; \db) = f^\circ(\xb^*; \db) + \mu h^\circ(\xb^*; \db) = M^\circ(\xb^*, \mu; \db) \geq 0 \, ,
\end{equation}

as desired.

§.§.§ Penalty approaches

In penalty based approaches, a direct-search scheme is applied to an auxiliary objective, that combines a measure of infeasibility with the objective function. A first approach of this form was the augmented Lagrangian technique
proposed by Lewis and Torczon [Lewis and Torczon, 2002], that handled nonlinear relaxable equality constraints as well as bounds on the objective. This technique was later extended to handle linear constraints in addition to bound
constraints [Kolda et al., 2006]. The underlying theoretical guarantees rely on using the stepsize as a stopping criterion for solving the subproblems.

In another line of work, Liuzzi et al. [Liuzzi and Lucidi, 2009] proposed an exact penalty approach to reformulate problem constrained with nonlinear, relaxable inequalities. The resulting penalized objective was defined on an
open set, and tackled with a line-search method akin to those described in Section <ref>. Liuzzi et al. [Liuzzi et al., 2010] considered a sequential penalty approach for problems with smooth objective,
inequality constraints and bound constraints. In this approach the inequality constraints other than bounds are relaxed and a line-search method is applied to minimize $f(\xb) + \frac{1}{\varepsilon} h_2(\xb)$ under bound constraints, where
$h_2$ is the smoothed infeasibility measure defined in (<ref>) and $\varepsilon>0$. A condition based on an upper bound on the stepsizes depending from $\mu$ is used to determine whether to decrease the penalty parameter $\varepsilon$
at each iteration. The case of nonsmooth relaxable constraints combined with bounds was addressed by Fasano et al. [Fasano et al., 2014]. Under constraint qualification and for small enough $\varepsilon$, it was then shown
that all the stationary points corresponding to minimizing $f(\xb) + \frac{1}{\varepsilon} h_1(\xb)$ under bound constraints are also (Clarke) stationary for the original problem.

§ CONCLUSION

Since their formal introduction more than 60 years ago, direct-search algorithms have been endowed with a rich convergence theory, that goes
significantly beyond a basic analysis for smooth, unconstrained problems. In addition to possessing guarantees in that setting, direct-search
methods can also be designed to cope with nonsmoothness in the objective function, constraints, and noisy estimates of all functions involved in
the problem formulation. In this survey, we provided key results that we hope will convince the readers of the strong theoretical foundations upon
which modern direct-search algorithms are built.

Our presentation has highlighted several results that we believe will play a key role in the future of direct-search techniques. Recent advances in
analyzing randomized direct-search techniques pave the way for scaling up those algorithms, possibly through combination with other tools from
high-dimensional statistics. Using modern constrained optimization techniques, such as Frank–Wolfe updates and manifold theory, provide
principled ways of handling common constraints and sparse structures. Stochastic and numerical noise in function values is likely to keep growing
in importance to address problems involving massive amounts of data samples such as hyperparameter tuning. These research directions are
likely to form the basis of next-generation direct-search methods.

[Achddou et al., 2022]
J. Achddou, O. Cappé, and A. Garivier.
Regret analysis of the stochastic direct search method for blind
resource allocation.
arXiv, 2022.
URL <https://arxiv.org/abs/2210.05222>.

[Alarie et al., 2021]
S. Alarie, C. Audet, P.-Y. Bouchet, and S. Le Digabel.
Optimisation of stochastic blackboxes with adaptive precision.
SIAM Journal on Optimization, 310 (4):0
3127–3156, 2021.

[Audet, 2014]
C. Audet.
A survey on direct search methods for blackbox optimization and their
In P. Pardalos and T. Rassias, editors, Mathematics Without
Boundaries: Surveys in Interdisciplinary Research, chapter 2, pages
31–56. Springer, New York, NY, 2014.

[Audet and Dennis, Jr., 2002]
C. Audet and J. E. Dennis, Jr.
Analysis of generalized pattern searches.
SIAM Journal on Optimization, 130 (3):0
889–903, 2002.

[Audet and Dennis, Jr., 2004]
C. Audet and J. E. Dennis, Jr.
A pattern search filter method for nonlinear programming without
SIAM Journal on Optimization, 140 (4):0
980–1010, 2004.

[Audet and Dennis, Jr., 2006]
C. Audet and J. E. Dennis, Jr.
Mesh adaptive direct search algorithms for constrained optimization.
SIAM Journal on Optimization, 170 (1):0
188–217, 2006.

[Audet and Dennis, Jr., 2009]
C. Audet and J. E. Dennis, Jr.
A progressive barrier for derivative-free nonlinear programming.
SIAM Journal on Optimization, 200 (1):0
445–472, 2009.

[Audet and Hare, 2017]
C. Audet and W. Hare.
Derivative-Free and Blackbox Optimization.
Springer Series in Operations Research and Financial Engineering.
Springer, Cham, Switzerland, 2017.

[Audet et al., 2010]
C. Audet, J. E. Dennis, Jr., and S. Le Digabel.
Globalization strategies for mesh adaptive direct search.
Computational Optimization and Applications, 460
(2):0 193–215, 2010.

[Audet et al., 2021]
C. Audet, K. J. Dzahini, M. Kokkolaras, and S. Le Digabel.
Stochastic mesh adaptive direct search for blackbox optimization
using probabilistic estimates.
Computational Optimization and Applications, 790
(1):0 1–34, 2021.

[Bäck, 1996]
T. Bäck.
Evolutionary Algorithms in Theory and Practice: Evolution
Strategies, Evolutionary Programming, Genetic Algorithms.
Oxford University Press, New York, NY, USA, 1996.

[Berahas et al., 2019]
A. S. Berahas, R. H. Byrd, and J. Nocedal.
Derivative-free optimization of noisy functions via quasi-newton
SIAM Journal on Optimization, 290 (2):0
965–993, 2019.

[Berahas et al., 2021]
A. S. Berahas, L. Cao, and K. Scheinberg.
Global convergence rate analysis of a generic line search algorithm
with noise.
SIAM Journal on Optimization, 310 (2):0
1489–1518, 2021.

[Bergou et al., 2020]
E. H. Bergou, E. Gorbunov, and P. Richtárik.
Stochastic three points method for unconstrained smooth minimization.
SIAM Journal on Optimization, 300 (4):0
2726–2749, 2020.

[Blanchet et al., 2019]
J. Blanchet, C. Cartis, M. Menickelly, and K. Scheinberg.
Convergence rate analysis of a stochastic trust-region method via
INFORMS Journal on Optimization, 10 (2):0
92–119, 2019.

[Bollapragada and Wild, 2023]
R. Bollapragada and S. M. Wild.
Adaptive sampling quasi-newton methods for zeroth-order stochastic
Mathematical Programming Computation, 150
(2):0 327–364, 2023.

[Bomze et al., 2021]
I. M. Bomze, F. Rinaldi, and D. Zeffiro.
Frank–wolfe and friends: a journey into projection-free first-order
optimization methods.
4OR, 190 (3):0 313–345, 2021.

[Booker et al., 1999]
A. J. Booker, J. E. Dennis, Jr., P. D. Frank, D. B. Serafini, V. Torczon, and
M. W. Trosset.
A rigorous framework for optimization of expensive functions by
Structural and Multidisciplinary Optimization, 170
(1):0 1–13, 1999.

[Bouchet, 2023]
P.-Y. Bouchet.
Théorie de l'optimisation sans dérivées dans le cas
PhD thesis, Polytechnique Montréal, December 2023.

[Boumal, 2023]
N. Boumal.
An introduction to optimization on smooth manifolds.
Cambridge University Press, Cambridge, United Kingdom, 2023.

[Brilli et al., 2023]
A. Brilli, M. Kimiaiei, G. Liuzzi, and S. Lucidi.
Worst case complexity bounds for linesearch-type derivative-free
arXiv, 2023.
URL <https://arxiv.org/abs/2302.05274>.

[Bueno et al., 2013]
L. F. Bueno, A. Friedlander, J. M. Martinez, and F. N. C. Sobral.
Inexact restoration method for derivative-free optimization with
smooth constraints.
SIAM Journal on Optimization, 230 (2):0
1189–1213, 2013.

[Cartis et al., 2022]
C. Cartis, N. I. M. Gould, and P. L. Toint.
Evaluation Complexity of Algorithms for Nonconvex Optimization:
Theory, Computation and Perspectives, volume MO30 of MOS-SIAM Series
on Optimization.
SIAM, 2022.

[Chang, 2012]
K. H. Chang.
Stochastic Nelder–Mead simplex method – A new globally
convergent direct search method for simulation optimization.
European Journal of Operational Research, 2200
(3):0 684–694, 2012.

[Chen et al., 2018]
R. Chen, M. Menickelly, and K. Scheinberg.
Stochastic optimization using a trust-region method and random
Mathematical Programming, 1690 (2):0 447–487,

[Clarke, 1983]
F. H. Clarke.
Optimization and Nonsmooth Analysis.
John Wiley and Sons, New York, 1983.
ISBN 9780898712568.
Reissued in 1990 by SIAM Publications, Philadelphia, as Vol. 5 in the
series Classics in Applied Mathematics.

[Conn et al., 2009]
A. R. Conn, K. Scheinberg, and L. N. Vicente.
Introduction to Derivative-Free Optimization.
SIAM, Philadelphia, 2009.
ISBN 978-0-898716-68-9.

[Coope and Price, 2001]
I. D. Coope and C. J. Price.
On the convergence of grid-based methods for unconstrained
SIAM Journal on Optimization, 110 (4):0
859–869, 2001.

[Cristofari and Rinaldi, 2021]
A. Cristofari and F. Rinaldi.
A derivative-free method for structured optimization problems.
SIAM Journal on Optimization, 310 (2):0
1079–1107, 2021.

[Custódio and Vicente, 2012]
A. L. Custódio and L. N. Vicente.
Analysis of direct searches for discontinuous functions.
Mathematical Programming, 133:0 299–325, 2012.

[De Leone et al., 1984]
R. De Leone, M. Gaudioso, and L. Grippo.
Stopping criteria for linesearch methods without derivatives.
Mathematical Programming, 30:0 285–300, 1984.

[Dennis, Jr. et al., 2004]
J. E. Dennis, Jr., C. J. Price, and I. D. Coope.
Direct search methods for nonlinearly constrained optimization using
filters and frames.
Optimization and Engineering, 50 (2):0
123–144, 2004.

[Diniz-Ehrhardt et al., 2008]
M. A. Diniz-Ehrhardt, J. M. Martínez, and M. Raydán.
A derivative-free nonmonotone line-search technique for unconstrained
Journal of computational and applied mathematics, 2190
(2):0 383–397, 2008.

[Diouane et al., 2023]
Y. Diouane, V. Kungurtsev, F. Rinaldi, and D. Zeffiro.
Inexact direct-search methods for bilevel optimization problems.
arXiv, 2023.
URL <https://arxiv.org/abs/2307.09924>.

[Dodangeh and Vicente, 2016]
M. Dodangeh and L. N. Vicente.
Worst case complexity of direct search under convexity.
Mathematical Programming, 1550 (1-2):0
307–332, 2016.

[Dodangeh et al., 2016]
M. Dodangeh, L. N. Vicente, and Z. Zhang.
On the optimal order of worst case complexity of direct search.
Optimization Letters, 10:0 699–708, 2016.

[Dolan et al., 2003]
E. D. Dolan, R. M. Lewis, and V. Torczon.
On the local convergence of pattern search.
SIAM Journal on Optimization, 140 (2):0
567–583, 2003.

[Dreisigmeyer, 2006]
D. W. Dreisigmeyer.
Equality constraints, riemannian manifolds and direct search methods.
Technical Report LA-UR-06-7406, Los Alamos National Laboratory, 2006.

[Dreisigmeyer, 2018]
D. W. Dreisigmeyer.
Direct search methods on reductive homogeneous spaces.
Journal of Optimization Theory and Applications, 1760
(3):0 585–604, 2018.

[Dzahini, 2022]
K. J. Dzahini.
Expected complexity analysis of stochastic direct-search.
Computational Optimization and Applications, 810
(1):0 179–200, 2022.

[Dzahini and Wild, 2022]
K. J. Dzahini and S. M. Wild.
Sparse Johnson–Lindenstrauss transforms and analysis of their
extreme singular values.
arXiv, 2022.
URL <https://arxiv.org/abs/2212.14858>.

[Dzahini et al., 2023]
K. J. Dzahini, M. Kokkolaras, and S. Le Digabel.
Constrained stochastic blackbox optimization using a progressive
barrier and probabilistic estimates.
Mathematical Programming, 1980 (1):0 675–732,

[Fasano et al., 2014]
G. Fasano, G. Liuzzi, S. Lucidi, and F. Rinaldi.
A linesearch-based derivative-free approach for nonsmooth constrained
SIAM Journal on Optimization, 240 (3):0
959–992, 2014.

[García-Palomares and Rodríguez, 2002]
U. M. García-Palomares and J. F. Rodríguez.
New sequential and parallel derivative-free algorithms for
unconstrained minimization.
SIAM Journal on Optimization, 130 (1):0
79–96, 2002.

[Garmanjani and Vicente, 2013]
R. Garmanjani and L. N. Vicente.
Smoothing and worst-case complexity for direct-search methods in
nonsmooth optimization.
IMA Journal of Numerical Analysis, 330 (3):0
1008–1028, 2013.

[Gasparo et al., 2001]
M. G. Gasparo, A. Papini, and A. Pasquali.
Nonmonotone algorithms for pattern search methods.
Numerical Algorithms, 28:0 171–186, 2001.

[Ghadimi and Lan, 2013]
S. Ghadimi and G. Lan.
Stochastic first-and zeroth-order methods for nonconvex stochastic
SIAM Journal on Optimization, 230 (4):0
2341–2368, 2013.

[Goldstein, 1977]
A. A. Goldstein.
Optimization for Lipschitz continuous functions.
Mathematical Programming, 130 (1):0 14–22,

[Gratton and Vicente, 2014]
S. Gratton and L. N. Vicente.
A merit function approach for direct search.
SIAM Journal on Optimization, 240 (4):0
1980–1998, 2014.

[Gratton et al., 2015]
S. Gratton, C. W. Royer, L. N. Vicente, and Z. Zhang.
Direct search based on probabilistic descent.
SIAM Journal on Optimization, 250 (3):0
1515–1541, 2015.

[Gratton et al., 2019]
S. Gratton, C. W. Royer, L. N. Vicente, and Z. Zhang.
Direct search based on probabilistic feasible descent for bound and
linearly constrained problems.
Computational Optimization and Applications, 720
(3):0 525–559, 2019.

[Grippo and Rinaldi, 2015]
L. Grippo and F. Rinaldi.
A class of derivative-free nonmonotone optimization algorithms
employing coordinate rotations and gradient approximations.
Computational Optimization and Applications, 60:0
1–33, 2015.

[Grippo et al., 1988]
L. Grippo, F. Lampariello, and S. Lucidi.
Global convergence and stabilization of unconstrained minimization
Journal of Optimization Theory and Applications, 56:0
385–406, 1988.

[Hooke and Jeeves, 1961]
R. Hooke and T. A. Jeeves.
“Direct search" solution of numerical and statistical problems.
Journal of the Association for Computing Machinery, 80
(2):0 212–229, 1961.

[Jin et al., 2021]
B. Jin, K. Scheinberg, and M. Xie.
High probability complexity bounds for line search based on
stochastic oracles.
In NeurIPS, volume 34, pages 9193–9203, 2021.

[Jin et al., 2023]
B. Jin, K. Scheinberg, and M. Xie.
Sample complexity analysis for adaptive optimization algorithms with
stochastic oracles.
arXiv, 2023.
URL <https://arxiv.org/abs/2303.06838>.

[Kelley, 2011]
C. T. Kelley.
Implicit Filtering.
Software Environment and Tools. Society for Industrial and Applied
Mathematics, Philadelphia, 2011.

[Kolda et al., 2003]
T. G. Kolda, R. M. Lewis, and V. Torczon.
Optimization by direct search: New perspectives on some classical and
modern methods.
SIAM Review, 450 (3):0 385–482, 2003.

[Kolda et al., 2006]
T. G. Kolda, R. M. Lewis, and V. Torczon.
Stationarity results for generating set search for linearly
constrained optimization.
SIAM Journal on Optimization, 170 (4):0
943–968, 2006.

[Kolda et al., 2006]
T. G. Kolda, R. M. Lewis, and V. Torczon.
A generating set direct search augmented lagrangian algorithm for
optimization with a combination of general and linear constraints.
Technical report, Technical Report SAND2006-5315, Sandia National
Laboratories, 2006.

[Konečnỳ and Richtárik, 2014]
J. Konečnỳ and P. Richtárik.
Simple complexity analysis of simplified direct search.
arXiv, 2014.
URL <https://arxiv.org/abs/1410.0390>.

[Kungurtsev et al., 2022]
V. Kungurtsev, F. Rinaldi, and D. Zeffiro.
Retraction based direct search methods for derivative free riemannian
arXiv, 2022.
URL <https://arxiv.org/abs/2202.11052>.

[Lacoste-Julien and Jaggi, 2015]
S. Lacoste-Julien and M. Jaggi.
On the global linear convergence of Frank–Wolfe optimization
Advances in Neural Information Processing Systems, 28, 2015.

[Larson et al., 2019]
J. Larson, M. Menickelly, and S. M. Wild.
Derivative-free optimization methods.
Acta Numerica, 28:0 287–404, 2019.

[Le Digabel, 2011]
S. Le Digabel.
Algorithm 909: NOMAD: Nonlinear optimization with the MADS
ACM Transactions on Mathematical Software, 370
(4):0 44:1–44:15, 2011.

[Le Digabel and Wild, 2023]
S. Le Digabel and S. M. Wild.
A taxonomy of constraints in black-box simulation-based optimization.
Optimization and Engineering, pages 1–19, 2023.

[Lewis and Torczon, 2002]
R. M. Lewis and V. Torczon.
A globally convergent augmented lagrangian pattern search algorithm
for optimization with general constraints and simple bounds.
SIAM Journal on Optimization, 120 (4):0
1075–1089, 2002.

[Lewis et al., 2000]
R. M. Lewis, V. Torczon, and M. W. Trosset.
Direct search methods: Then and now.
Journal of Computational and Applied Mathematics, 1240
(1–2):0 191–207, 2000.

[Lin et al., 2022]
T. Lin, Z. Zheng, and M. Jordan.
Gradient-free methods for deterministic and stochastic nonsmooth
nonconvex optimization.
Advances in Neural Information Processing Systems,
35:0 26160–26175, 2022.

[Liuzzi and Lucidi, 2009]
G. Liuzzi and S. Lucidi.
A derivative-free algorithm for inequality constrained nonlinear
programming via smoothing of an $\ell_\infty$ penalty function.
SIAM Journal on Optimization, 200 (1):0 1–29,

[Liuzzi et al., 2010]
G. Liuzzi, S. Lucidi, and M. Sciandrone.
Sequential penalty derivative-free methods for nonlinear constrained
SIAM Journal on Optimization, 200 (5):0
2614–2635, 2010.

[Lucidi and Sciandrone, 2002]
S. Lucidi and M. Sciandrone.
A derivative-free algorithm for bound constrained optimization.
Computational Optimization and Applications, 210
(2):0 119–142, 2002a.

[Lucidi and Sciandrone, 2002]
S. Lucidi and M. Sciandrone.
On the global convergence of derivative-free methods for
unconstrained optimization.
SIAM Journal on Optimization, 130 (1):0
97–116, 2002b.

[Lucidi et al., 2002]
S. Lucidi, M. Sciandrone, and P. Tseng.
Objective-derivative-free methods for constrained optimization.
Mathematical Programming, 92:0 37–59, 2002.

[Nelder and Mead, 1965]
J. A. Nelder and R. Mead.
A simplex method for function minimization.
The computer journal, 70 (4):0 308–313, 1965.

[Nesterov and Spokoiny, 2017]
Y. Nesterov and V. Spokoiny.
Random gradient-free minimization of convex functions.
Foundations of Computational Mathematics, 17:0
527–566, 2017.

[Paquette and Scheinberg, 2020]
C. Paquette and K. Scheinberg.
A stochastic line search method with expected complexity analysis.
SIAM Journal on Optimization, 300 (1):0
349–376, 2020.

[Pauwels, 2016]
B. Pauwels.
Optimisation sans dérivées sous incertitudes appliquées
à des simulateurs coûteux.
PhD thesis, Université de Toulouse, March 2016.
URL <http://thesesups.ups-tlse.fr/3103/1/2016TOU30035.pdf>.

[Rinaldi and Zeffiro, 2023]
F. Rinaldi and D. Zeffiro.
Avoiding bad steps in Frank–Wolfe variants.
Computational Optimization and Applications, 840
(1):0 225–264, 2023.

[Rinaldi et al., 2023]
F. Rinaldi, L. N. Vicente, and D. Zeffiro.
Stochastic trust-region and direct-search methods: A weak tail bound
condition and reduced sample sizing.
arXiv, 2023.
URL <https://arxiv.org/abs/2202.11074>.

[Roberts and Royer, 2023]
L. Roberts and C. W. Royer.
Direct search based on probabilistic descent in reduced spaces.
SIAM Journal on Optimization, 330 (4):0
3057–3082, 2023.

[Torczon, 1997]
V. Torczon.
On the convergence of pattern search algorithms.
SIAM Journal on Optimization, 70 (1):0 1–25,

[Vicente, 2013]
L. N. Vicente.
Worst case complexity of direct search.
EURO Journal on Computational Optimization, 10
(1):0 143–153, 2013.

[Wright, 1996]
M. H. Wright.
Direct search methods: Once scorned, now respectable.
Pitman Research Notes in Mathematics Series, 344:0
191–208, 1996.
URL <https://api.semanticscholar.org/CorpusID:64988980>.

[Wright and Recht, 2022]
S. J. Wright and B. Recht.
Optimization for Data Analysis.
Cambridge University Press, 2022.

The submitted manuscript has been co-created by UChicago Argonne, LLC, Operator of
Argonne National Laboratory (“Argonne”). Argonne, a U.S. Department of
Energy Office of Science laboratory, is operated under Contract No.
The U.S. Government retains for itself, and others acting on its behalf, a
paid-up nonexclusive, irrevocable worldwide license in said article to
reproduce, prepare derivative works, distribute copies to the public, and
perform publicly and display publicly, by or on behalf of the Government.  The
Department of Energy will provide public access to these results of federally
sponsored research in accordance with the DOE Public Access Plan.


# HEROES: The Hawaii eROSITA Ecliptic Pole Survey Catalog

A. J. Taylor Department of Astronomy, University of Wisconsin–Madison, 475 N.
Charter Street, Madison, WI 53706, USA  A. J. Barger Department of Astronomy,
University of Wisconsin–Madison, 475 N. Charter Street, Madison, WI 53706, USA
Department of Physics and Astronomy, University of Hawaii, 2505 Correa Road,
Honolulu, HI 96822, USA  Institute for Astronomy, University of Hawaii, 2680
Woodlawn Drive, Honolulu, HI 96822, USA  L. L. Cowie Institute for Astronomy,
University of Hawaii, 2680 Woodlawn Drive, Honolulu, HI 96822, USA  G.
Hasinger European Space Astronomy Centre (ESA/ESAC), E-28691 Villanueva de la
Cañada, Madrid, Spain  E. M. Hu  Institute for Astronomy, University of
Hawaii, 2680 Woodlawn Drive, Honolulu, HI 96822, USA  A. Songaila  Institute
for Astronomy, University of Hawaii, 2680 Woodlawn Drive, Honolulu, HI 96822,
USA

###### Abstract

We present a seven band ($g$, $r$, $i$, $z$, $y$, NB816, NB921) catalog
derived from a Subaru Hyper Suprime-Cam (HSC) imaging survey of the North
Ecliptic Pole (NEP). The survey, known as HEROES, consists of 44 deg2 of
contiguous imaging reaching median 5$\sigma$ depths of $g$: 26.5, $r$: 26.2,
$i$: 25.7, $z$: 25.1, $y$: 23.9, NB816: 24.4, NB921: 24.4 mag. We reduced
these data with the HSC pipeline software hscPipe, and produced a resulting
multiband catalog containing over 25 million objects. We provide the catalog
in three formats: (1) a collection of hscPipe format forced photometry
catalogs, (2) a single combined catalog containing every object in that
dataset with selected useful columns, and (3) a smaller variation of the
combined catalog with only essential columns for basic analysis or low memory
machines. The catalog uses all the available HSC data on the NEP and may serve
as the primary optical catalog for current and future NEP deep fields from
instruments and observatories such as SCUBA-2, eROSITA, Spitzer, Euclid, and
JWST.

Catalogs (205), Galaxy counts (588), Surveys (1671), Observational cosmology
(1146)

††facilities: Subaru Telescope††software: astropy: Astropy Collaboration et
al. (2013, 2018), hscPipe: Bosch et al. (2018)

##  1 Introduction

Since its first light in 2013, Hyper Suprime-Cam (HSC; Miyazaki et al., 2018)
on the Subaru 8.2m Telescope has been the premier wide field optical imager on
6–10 meter class telescopes. The large collecting area of Subaru and the wide
field of view of HSC (1$\fdg$5) allow for the efficient observation of both
wide and deep surveys. The largest of these projects—the HSC Subaru Strategic
Program (HSC-SSP)—has fully mapped 670 deg2 of sky in $grizy$ broadband
filters and over 1470 deg2 in partially observed (incomplete coverage in
$grizy$) area (Aihara et al., 2018, 2019, 2022). HSC-SSP has also produced
deep observations in the DEEP2-3, SXDS+XMM-LSS, ELAIS-N1, and COSMOS fields
totaling $>$30 deg2 of imaging in $grizy$ broadband filters, as well as NB387,
NB816, and NB921 narrowband filters, all at a 5$\sigma$ depth $>$25 mag.

In 2015, given the enormous potential of HSC and its available narrowband
filters, we designed and began to execute HEROES: the Hawaii EROsita Ecliptic
pole Survey. Our goal was to produce a survey of the North Ecliptic Pole (NEP)
that would match the filter coverage and imaging depths of the HSC-SSP Deep
fields with a much wider contiguous area, as well as complement the then-
future deepest eROSITA X-ray observations (Merloni et al., 2012).

In recent years, the NEP has become a major focus of a number of ground and
space-based surveys and missions spanning sub-millimeter to X-ray wavelengths.
HEROES provides complementary wide-field optical broadband and narrowband data
to these surveys.

At sub-millimeter wavelengths, the NEP has also been extensively observed with
SCUBA-2 on the James Clerk Maxwell Telescope through the S2CLS (0.6 deg2, 850
$\mu$m; Geach et al., 2017), NEPSC2 (2 deg2, 850 $\mu$m; Shim et al., 2020),
and S2TDF (0.087 deg2, 850$\mu$m; Hyun et al., 2023) surveys. HEROES will
provide sub-millimeter studies with corresponding optical data that may be
used to find optical counterparts for sub-millimeter bright sources. The
combination of optical/NIR and sub-millimeter flux may better constrain
photometric redshifts and spectral energy distribution fittings for the
determination of stellar mass, star formation rate, age, and dust attenuation
for such sources (e.g. S. McKay et al. submitted).

In the infrared, the NEP contains the Spitzer IRAC Dark Field (Krick et al.,
2008, 2009; Frost et al., 2009), as well as the upcoming 20 deg2 Euclid Deep
Field North (Amendola et al., 2013, 2018). HEROES will serve as a natural
complement and extension of these datasets into optical wavelengths, providing
broadband coverage from 0.4–0.8 $\mu$m (when combined with Spitzer/IRAC) and
providing individual $griz$ magnitudes to help better complement the upcoming
single 550–900 nm very-broadband Euclid/VIS observations and 0.92–2 $\mu$m
Euclid/NISP $YJH$ photometric and 1.1–2 $\mu$m slitless spectroscopic
observations (Laureijs et al., 2011). The HEROES data will permit improved SED
fitting for Euclid detected sources and ultimately provide a more robust
understanding of the properties of these objects.

In x-rays, the NEP is already home to the deepest eROSITA X-ray observations
(Merloni et al., 2012), and will contain the future SPHEREx Deep North Field
(Doré et al., 2016, 2018). For these missions, HEROES and may provide insight
into optical properties of x-ray detected AGN including photometric redshifts,
as well as provide superior target astrometry when compared to the x-ray
measurements. For example, Radzom et al. (2022) used Chandra data in
conjunction with HSC data in the SSA22 field to produce x-ray luminosity
functions at redshifts $z=0.2--4$. With the combination of HEROES, eROSITA,
SPHEREx, and Euclid, this type of study could be replicated in the NEP with
$>500\times$ the survey area.

Moreover, HEROES and the NEP are especially important for space-based
observatories that orbit the Sun-Earth L2 Lagrange point, as the NEP and the
South Ecliptic Pole are typically part of any such observatory’s continuous
viewing zone. As such, the NEP contains the aforementioned Sptizer IRAC Dark
Field, eROSITA Deep Field, upcoming Euclid Deep Field North, and SPHEREx Deep
Field North as well as the JWST Time Domain Field (TDF; part of the PEARLS
project; Windhorst et al., 2017, 2022).

HEROES was initially reduced in 2017 using the Pan-STARRS Image Processing
Pipeline (IPP; Magnier et al., 2016, 2020a, 2020b). This version of HEROES had
incomplete coverage in the $r$ and $y$ bands at low R.A. pointings and did not
have the additional JWST TDF pointing (see §2 below).

HEROES is the largest contiguous HSC narrowband survey to date. As such, we
previously used the wide-field narrowband coverage in the initial HEROES
dataset for studies of Lyman-$\alpha$ Emitters (LAEs) near the epoch of
reionization at $z=5.7-6.6$ (Songaila et al., 2018, 2022; Taylor et al., 2020,
2021), as well as the development of a broadband selection technique for
emission line galaxies (Rosenwasser et al., 2022).

In 2021, we completed the final HEROES observations. Here we present the
photometric catalog of the complete and newly re-reduced version of HEROES,
which also incorporates all the archival HSC data on this field.

In Section 2, we describe the HSC data. In Section 3, we summarize the data
reduction and processing. In Section 4, we present the final catalog and
describe its format and availability. In Section 5, we verify the catalog’s
quality and measure its depth in each filter. Finally, in Section 6, we
demonstrate both a $z=6.6$ LAE sample selection and $g$, $r$, and $i$-band
dropout selections.

We assume $\Omega_{M}=0.3$, $\Omega_{\Lambda}=0.7$, and $H_{0}=70$ km s-1
Mpc-1 throughout. All magnitudes are given in the AB magnitude system, where
an AB magnitude is defined by $m_{AB}=-2.5\log f_{\nu}-48.60$. Here $f_{\nu}$
is the flux of the source in units of ergs cm-2 s-1 Hz-1.

##  2 Observations

We summarize the NEP HSC observations in Table 1 and illustrate the pointing
centers in Figure 1. Our HEROES observations contributed 2574 targeted shots
(a “shot” is a single HSC exposure consisting of 112 mosaicked CCDs). We used
a hexagonal grid of 38 pointings separated by 1$\fdg$0 to balance pointing
overlap with the 1$\fdg$5 diameter HSC field of view (see Figure 1). In our
pointings, we used a six point mosaic dither pattern with one central shot and
five shots evenly distributed around a 2$\farcm$0 radius circle centered on
the central shot (see Songaila et al., 2018, their Figure 1). In addition to
our pointing grid, we conducted observations at a specialized pointing with
the HSC-r2, HSC-z, and NB921 filters centered on the JWST TDF (Windhorst et
al., 2017, 2022). Our HEROES observations were taken using the standard set of
HSC broadband filters: HSC-g, HSC-r2, HSC-i2, HSC-z, HSC-Y, and the NB816 and
NB921 narrowband filters (hereafter, referred to as $g$, $r$, $i$, $z$, $y$,
NB816, and NB921).

In addition to the HEROES NEP observations, the AKARI-HSC survey (Oi et al.,
2021, hereafter, AKARI) previously observed 5.4 deg2 of the NEP in $grizy$
broadbands to complement the infrared coverage of the AKARI satellite
(Matsuhara et al., 2006; Murakami et al., 2007; Lee et al., 2009). The Hawaii
Two-0 (H20; Beck et al., 2020) survey also combined parts of the HEROES raw
data with new observations to provide $griz$ broadband coverage to complement
the upcoming Euclid Deep Field North (Amendola et al., 2013, 2018). We
incorporated 181 archival shots in HSC $g$, $z$, and $y$ filters from the
AKARI (Oi et al., 2021) and 305 archival shots in $g$, $r$, $i$, and $z$
filters from H20 (Beck et al., 2020) into our present HEROES data reduction.
All together, these observations cover 44 deg2 with full 7-band imaging.

Figure 1: Pointing centers and overlap for the HSC data in the NEP. Our HEROES pointings are shown in blue, the H20 pointings in red, and the AKARI pointings in green. Our pointing targeting the JWST TDF is shown in yellow. Each of the pointings used either the HEROES dither patten or something similar (see §2). Color shading shows pointing overlaps but does not explicitly indicate stacked imaging depth. Table 1: HSC Observations FilteraaWe refer to the HSC-g, HSC-r2, HSC-i2, HSC-z, and HSC-Y filters as $g$, $r$, $i$, $z$, and $y$ respectively throughout this article. | HEROES | H20 | AKARI | Total
---|---|---|---|---
| Shots | Hours | Shots | Hours | Shots | Hours | Shots | Hours
HSC-g | 355 | 8.90 | 40 | 2.78 | 100 | 7.92 | 495 | 19.59
HSC-r2 | 208 | 9.02 | 60 | 3.33 |  $\cdots$ |  $\cdots$ | 268 | 12.35
HSC-i2 | 452 | 16.80 | 100 | 6.14 |  $\cdots$ |  $\cdots$ | 552 | 22.94
HSC-z | 735 | 24.21 | 105 | 7.06 | 29 | 1.85 | 869 | 33.11
HSC-Y | 241 | 10.78 |  $\cdots$ |  $\cdots$ | 52 | 4.55 | 293 | 12.23
NB816 | 235 | 11.33 |  $\cdots$ |  $\cdots$ |  $\cdots$ |  $\cdots$ | 235 | 11.33
NB921 | 348 | 20.86 |  $\cdots$ |  $\cdots$ |  $\cdots$ |  $\cdots$ | 348 | 20.86
Total | 2574 | 101.79 | 305 | 19.31 | 181 | 14.32 | 3060 | 135.42

##  3 Data Reduction

We processed the HEROES dataset using the standard procedure provided by the
hscPipe pipeline version 8.4 (Bosch et al., 2018). We performed this
processing on the National Astronomical Observatory of Japan Large-scale Data
Analysis System computing cluster with access provided from our September 2021
HSC observations. hscPipe is a comprehensive pipeline for the reduction and
processing of HSC data that acts in 5 main steps: Bias/Dark/Flat/Fringe
processing, Single visit processing, Moasicking, Coadding images, and
Multiband analysis. In each of these steps, we used the default hscpipe
parameters, unless noted below.

We assembled preprocessed Bias/Dark/Flat/Fringe data from the HSC calibration
data archive111https://www.naoj.org/Observing/Instruments/HSC/calib_data.html
\- last accessed 2023 February 13 that matched all of the science data filters
and observing runs. We performed single visit processing (detrending, WCS, and
photometric calibration) on all individual CCD frames (103 science frames per
shot) using the pipeline defaults. In the mosaicking step, we defined a custom
tract and patch scheme. A hscPipe tract is an area of sky over which images
are mosaicked and coadded with a common WCS solution. To produce data products
of manageable size, a tract is divided into a grid of images called patches.
We defined a single tract that encompassed the entire data coverage, and we
divided it into a grid of $20\times 14$ patches, each measuring 10,200
$\times$ 10,200 pixels to allow for overlap between adjacent patches. We show
the position of the patches on the sky in Figure 2. We used the hscPipe FGCM
(Forward Global Calibration Method) to calibrate the data simultaneously in
all seven filters to Pan-STARRS photometry and astrometry.

Figure 2: Patch locations for the HSC data in the NEP. The bounding box of
each patch is shown in black and labeled with its patch ID. The field imaging
coverage is shown in the background in light gray shading.

We executed the standard imaging coadding routines for each filter. We then
ran the standard multiband analysis pipeline procedures to combine the
detected source catalogs across the seven different filters for each patch and
to perform forced photometry on the resulting combined catalogs. After
completing this step, we noticed that the coadded images at the edges of the
field included extrapolated background modeled pixels that extended beyond the
science imaging coverage. These extrapolated pixels resulted in erroneous
measurements of object fluxes during the multiband analysis step on our
initial processing run. To remedy this, we masked all pixels outside of the
science imaging coverage. We re-ran the multiband measurements after this
manual change and confirmed that our modification resolved the irregular
photometric measurements.

##  4 The Catalog

We present the final HEROES catalog in three forms to best provide the maximum
amount of available information, as well as convenient and manageable file
sizes.

First, we provide the forced photometry output catalogs from the multiband
analysis step of hscPipe. These are organized by filter and patch using the
naming format forced_src-{filter}-0-{patch}.fits. The columns of these FITS
tables are provided in the FITS headers and are summarized in the hscPipe
documentation222https://hsc.mtk.nao.ac.jp/pipedoc/pipedoc_8_e/tutorial_e/schema_multiband.html
\- last accessed 2023 February 13. These 1671 catalogs each contain 239
columns and total 211 GB (compressed). While these catalogs provide all of the
available hscPipe forced photometry information, they are sometimes
inconvenient for use in multifilter full field studies and contain many
columns that are not useful for typical research applications.

We also provide two catalogs that each include all of the objects in the
dataset and select columns from the forced_src catalogs. The first catalog
(HEROES_Full_Catalog.fits, 25,445,387 objects, 172 columns, 12.4 GB) contains
all of the columns described in the Appendix in Table 2. The second catalog
(HEROES_Small_Catalog.fits, 25,445,387 objects, 37 columns, 4.1 GB) is
designed for more basic analyses or for machines with less system memory and
only contains a subset of the full catalog’s columns. It contains selected
columns, as noted in the Appendix in Table 2.

In these catalogs, we converted forced_src fluxes to AB magnitudes. To
preserve the information for sources with negative measured fluxes, we report
these values as negative magnitudes. For example; we report an object with a
measured flux of $-3631\times 10^{-11}$ Jy as an AB magnitude of -27.5. For
objects that do not have imaging coverage in a given filter, or lack measured
fluxes in a given filter, we report magnitudes of -99.

These data products are all publicly available at Harvard
Dataverse333https://dataverse.harvard.edu/dataverse/heroes.

##  5 Data Quality

We tested the processed data quality by calculating 2$\farcs$0 diameter
aperture magnitude depths across the field and in each filter with two
different methods. In our first method, we placed 100,000 apertures randomly
in each patch ($\sim$150 apertures per arcmin2) and measured the flux in each
aperture. We then divided the patch into 100 subregions and analyzed the
apertures in each region separately. We used sigma-clipping to discard
apertures that captured flux from field objects and took the standard
deviation of the flux measurements from the remaining apertures. We converted
this flux standard deviation to a magnitude to produce a 1$\sigma$ aperture
magnitude depth, which we then converted to a 5$\sigma$ depth for the
subregion by subtracting 1.75 mag ($-2.5\log_{10}5$). We report the overall
median 2$\farcs$0 diameter aperture magnitude 5$\sigma$ depths as: $g$: 26.0,
$r$: 25.7, $i$: 25.2, $z$: 24.7, $y$: 23.7, NB816: 24.3, NB921: 24.2.

In our second calculation method, we adopted the methodology of Oi et al.
(2021). Here, we divided the survey into $150\times 150$ subregions. In each
subregion, we used the catalog provided fluxes and flux errors to select
objects with a signal-to-noise of $\sim$5; ($5\pm 0.1$). We then took the
median magnitude of this S/N$\sim$5 population as the 5$\sigma$ depth for the
subregion. For this method, we report the overall median 2$\farcs$0 diameter
aperture magnitude 5$\sigma$ depths as: $g$: 26.5, $r$: 26.2, $i$: 25.7, $z$:
25.1, $y$: 23.9, NB816: 24.4, NB921: 24.4, and we show their variations across
the field in Figure 3. We attribute the 0.1–0.5 magnitude differences between
the two methods to potential pattern noise effects or correlated variance that
may not be fully sampled by the measurement routines in hscPipe, and/or to not
fully cleaning contamination in the aperture method. The two methods may
bracket the true limits.

Figure 3: 5$\sigma$ detection limits for each filter across the survey field
from the Oi et al. (2021) method. Note the depth increases in the centers of
the $g,r,i,z,y$ imaging from the overlapping archival data and the significant
depth increases in $r$, $z$, and NB921 near the JWST TDF (R.A. 17:22:47.896,
Decl. +65:49:21.54) from our additional targeted pointing.

We tested the data quality by comparing the spatial number density of sources
as a function of magnitude against other large surveys. To limit contaminating
sources in our catalog for these comparisons, we applied a number of cuts to
the catalog. First, we required all catalog objects (both stars and galaxies)
to have the is_primary flag. This flag selects objects that are not detected
as blended composite objects. For each filter, we also removed all objects
with the flag {filter}_base_PixelFlags_flag_edge. This cut has two main
effects: First, it removes all sources that are outside of the science imaging
coverage. Second, it removes bad sources near saturated pixels (e.g., near
diffraction spikes from bright stars). After these cuts, we compared the area
density of sources in our catalog to the HSC COSMOS Deep/UltraDeep (D/UD)
catalog from HSC-SSP (10.0 deg2, Aihara et al., 2021) and the Dark Energy
Survey Year 3 GOLD catalog (DESY3; 5347 deg2; Sevilla-Noarbe et al., 2021;
Hartley et al., 2022). We use 2$\farcs$0 diameter aperture magnitudes from
hscPipe for HEROES and COSMOS, and we use “Single Object Fitted Corrected”
Magnitudes for DESY3 (2$\farcs$0 diameter aperture magnitudes are not provided
in the DESY3 Gold Catalog). We show this comparison in Figure 4.

Figure 4: Area density of sources in the $g$, $r$, $i$, and $z$ bands as a
function of magnitude from HEROES (blue, 2” diameter aperture magnitudes),
COSMOS D/UD (red, 2” diameter aperture magnitudes, Aihara et al., 2021), and
the DESY3 GOLD catalog (green, single object fitted magnitudes, Sevilla-Noarbe
et al., 2021; Hartley et al., 2022). While the DESY3 catalog appears to be
slightly overdense in the $g$ and $r$ bands, HEROES shows strong agreement
with the area densities from COSMOS D/UD in all four bands and DESY3 in $i$
and $z$.

Across $g$, $r$, $i$, and $z$, HEROES (blue) shows excellent agreement with
COSMOS D/UD (red). The DESY3 Gold catalog appears to be slightly overdense in
the bluer $g$ and $r$ bands when compared to the other catalogs, but it shows
good agreement in the $i$ and $z$ bands. Differences between the catalogs at
the bright end are likely due to minor contamination by bright stars and image
artifacts, while differences between the catalogs at the faint end are due to
the different depths between the surveys. From these comparisons, we conclude
that HEROES is consistent with other leading surveys that have well-calibrated
photometry and low contamination.

We further characterize the contamination rate in HEROES through visual
inspections of 1000 sources drawn randomly from the filtered sample described
above. In these inspections, we inspect $grizy$ thumbnails of the sources and
look for diffraction spikes, glints, halos, or other visual artifacts that are
detected as the source under inspection. Of the 1000 inspected sources, we
find that only 27 (2.7$\pm$0.5%) are impacted by visual artifacts. In most
cases, the artifacts are the unsaturated tails of diffraction spikes or halo-
like arcs from bright stars in the field.

##  6 Sample Selections

Here we demonstrate as examples two different sample selections using HEROES.
These are the focus of our previous and upcoming research projects.

###  6.1 Narrowband Selection of $z\sim 6$ LAEs

We used a previous reduction of the HEROES data for our sample selections of
$z\sim 5$ and $z\sim 6$ LAEs in Taylor et al. (2020, 2021) and Songaila et al.
(2022). We now repeat these selections to verify the photometric consistency
between the two reductions.

Our selection criteria are as follows. First, we remove sources with bad
pixel, bright object, edge pixel, saturated center, cosmic ray, and
interpolated center flags in the $i$, $z$, $y$, and NB921 filters. We then
remove sources within 6$\farcs$0 of any Gaia star with Gaia $g$ magnitude less
than 8. In Figure 5, we show $1\%$ of the resulting clean sample as black
points. From this clean sample, we next require strong detections (n921_detect
= True, n921_base_SdssCentroid_flag = False, n921_base_SdssShape_flag = False)
in NB921 and 2$\sigma$ non-detections ({filter}_detect = False) with forced
aperture magnitude uncertainties greater than 0.5 mags in the $g$, $r$, and
$i$ bands to enforce a strong Lyman break blueward of the Ly$\alpha$ line. We
select a narrowband excess $z$ – NB921 $>1.3$ mags. We also adopt the $\Sigma$
parameter from Sobral et al. (2013). This parameter characterizes the
significance of a narrowband excess above the uncertainties in the NB921 and
$z$ source magnitudes and is given by

$\Sigma=\frac{1-10^{-0.4(z-NB921)}}{10^{-0.4(27-NB921)}\sqrt{\sigma_{NB921}^{2}+\sigma_{z}^{2}}}\,,$
(1)

where $z$ and NB921 are the AB magnitudes of $z$ and NB921, $\sigma_{\rm
NB921}$ and $\sigma_{z}$ are the average $1\sigma$ image count rate
uncertainties in 2$\farcs$0 diameter apertures in NB921 and $z$, and 27 is the
magnitude zeropoint of the imaging data. In our source selection, we require
$\Sigma>3$. We show the resulting cut in Figure 5 (blue curve).

We then visually inspect cutouts of the remaining sources in stacked $gri$,
$z$, $y$, and NB921 to reject sources with significant contamination from an
elevated background, glints, diffraction spikes, transients, or other
artifacts. This visual inspection is also helpful in rejecting sources that
are not detected in $g$, $r$, or $i$ separately but are visually identifiable
in stacked $gri$ cutouts. For a narrowband selection targeting objects with
$z$ – NB921 $>1.3$ and NB921$<24.25$, the above cuts produced a sample of 384
LAE candidates. After visual inspection, we reduced this sample to 63
candidates that showed no hint of emission in stacked and smoothed $gri$
cutouts, had compact morphologies, and were visually free of contamination.
This significant reduction in candidates through visual inspection is
primarily due to the ability of stacked $gri$ cutouts to detect low-redshift
sources that do not show significant $gri$ emission in single filter
observations. Furthermore, narrowband excess and Lyman break samples are more
susceptible to objects with visual artifacts and contamination, as many forms
of contamination may artificially emulate the narrowband excess criteria and
non-detections in the bluer bands. Removing these sources through stricter
magnitude and color cuts may risk reducing the sample completeness of the
inherently rare bright $z>6$ LAEs, thus we use visual inspections to eliminate
contaminating objects to ensure that our samples remain both complete and
pure.

From these selection criteria (excluding the sources at NB921$<24.25$ that
were not visually re-inspected), we completely recover the Taylor et al.
(2020) and Songaila et al. (2022) samples (shown in Figure 5 as red crosses)
and identify additional candidates for spectroscopic follow up (A. Songaila et
al. 2023, in prep). These candidates and the recovered previous samples are
uniformly distributed across the survey field with no obvious visual
clustering or gradient beyond minor correlations with the imaging depth.

Figure 5: Color excess vs. magnitude diagram for $z\sim 6.6$ LAEs in HEROES.
The black points are 1% of the HEROES catalog that pass our initial source
quality cuts. The blue curve is the narrowband excess significance cut (Sobral
et al., 2013), and the dashed black line is the narrowband excess cut. The red
crosses are the spectroscopically confirmed $z\sim 6.6$ LAEs from Taylor et
al. (2020) and Songaila et al. (2022).

###  6.2 Broadband Selection of Dropout Galaxies

In order to test further the quality and science potential of the dataset, we
also demonstrate a broadband dropout selection using the selection criteria
and color-color cuts from Ono et al. (2018). In their study, “GOLDRUSH”, they
selected $z\sim 4,5,6,7$ galaxies using the dropout method with the UltraDeep,
Deep, and Wide fields from HSC-SSP. These fields total 102.7 deg2 in combined
area, with the largest field (W-XMM) providing 28.5 deg2 of coverage. We are
currently working on a full comparison with the GOLDRUSH luminosity functions
and clustering analysis (Harikane et al., 2018), and we summarize the
preliminary galaxy selection results below.

As both the HEROES and HSC-SSP catalogs are produced by hscpipe, it is simple
to adopt the Ono et al. (2018) selection criteria from their Table 2. In
brief, they first required sources to have no bad pixel, bright object, edge
pixel, saturated center, cosmic ray, or interpolated center flags in $grizy$.
For each sample of $g$, $r$, $i$, and $z$-dropouts, they required non-
detections in filters blueward of the dropout filter and strong detections in
filters redward of the dropout filter. They then adopted color-color cuts from
Hildebrandt et al. (2009) (see Ono et al., 2018, Equations 1–10) to produce
their initial dropout selections. We adopt the same cuts and show the color-
color criteria and our results in Figure 6.

Figure 6: Left: 1% of the objects that passed the quality cuts but not the
$g$-dropout color-color cuts (black points), and 1% of the objects that also
passed the $g$-dropout color-color cuts (red points). Center: 1% of the
objects that passed the quality cuts but not the $r$-dropout color-color cuts
(black points), and 10% of the objects that also passed the $r$-dropout color-
color cuts (red points). Right: Center: 1% of the objects that passed the
quality cuts but not the $i$-dropout color-color cuts (black points), and 100%
of the objects that also passed the $i$-dropout color-color cuts (red points).

In each panel of Figure 6, we show a subset of the sources that pass the above
described quality cuts as black points, and those that also pass the color-
color dropout cuts as red points. For $g$, $r$, and $i$ dropouts, we find
295129, 18607, and 124 galaxies, respectively. This corresponds to $g$: 6700,
$r$: 420, and $i$: 2.8 dropouts deg-2. These surface densities are comparable
to the densities from Ono et al. (2018) of $g$: 5300, $r$: 380, $i$: 5.2
deg-2, and we attribute the offsets to differences in imaging depth and
Poisson statistics.

The surface densities of all three classes of dropouts are roughly uniform,
with differences of no more than a factor of 2 over the HEROES field due
primarily to differences in imaging depth both between bands and across the
field. We will refine these selections and compare our resulting luminosity
functions and galaxy-galaxy angular correlation functions to the GOLDRUSH
sample in A. J. Taylor et al. 2023, in prep.

##  7 Summary

We present the complete photometric catalog from HEROES: a 44 deg2 Subaru/HSC
imaging survey of the NEP in $grizy$ broadbands and NB816+NB921 narrowbands.
The catalog contains 25.4 million objects and is available in patch by patch,
filter by filter hscpipe forced_src format, as well as in two combined
catalogs with selected columns.

HEROES has enormous potential due to its overlap with other legacy, current,
and future missions and surveys (e.g., AKARI, eROSITA, H20, S2CLS, NEPSC2,
Spitzer, Euclid, JWST TDF). Outside of these complementary datasets, we are
using HEROES to produce luminosity functions and angular correlation functions
for $z\sim 3-7$ dropout galaxies, as well as continuing to search for LAEs
near the epoch of reionization. We hope this public catalog release will
enable new studies of galaxy evolution across cosmic time and provide
complementary optical data for upcoming NEP surveys.

We thank the anonymous referee for their constructive report that helped us to
improve this work. We gratefully acknowledge support for this research from
NSF grants AST-1715145 (A. J. B) and AST-1716093 (E. M. H., A. S.). We also
gratefully acknowledge the William F. Vilas Estate (A. J. T.) and a Kellett
Mid-Career Award and a WARF Named Professorship from the University of
Wisconsin-Madison Office of the Vice Chancellor for Research and Graduate
Education with funding from the Wisconsin Alumni Research Foundation (A. J.
B.). This research is based on data collected at the Subaru Telescope, which
is operated by the National Astronomical Observatory of Japan. Data analysis
was carried out on the Multi-wavelength Data Analysis System operated by the
Astronomy Data Center and the Large-scale Data Analysis System co-operated by
the Astronomy Data Center and the Subaru Telescope. We especially thank the
HSC Software Help Desk for their rapid replies and helpful support during the
data reduction process. We wish to recognize and acknowledge the very
significant cultural role and reverence that the summit of Maunakea has always
had within the indigenous Hawaiian community. We are most fortunate to have
the opportunity to conduct observations from this mountain.

## References

  * Aihara et al. (2018) Aihara, H., Armstrong, R., Bickerton, S., et al. 2018, PASJ, 70, S8
  * Aihara et al. (2019) Aihara, H., AlSayyad, Y., Ando, M., et al. 2019, PASJ, 71, 114
  * Aihara et al. (2021) —. 2021, arXiv e-prints, arXiv:2108.13045
  * Aihara et al. (2022) —. 2022, PASJ, 74, 247
  * Amendola et al. (2013) Amendola, L., Appleby, S., Bacon, D., et al. 2013, Living Reviews in Relativity, 16, 6
  * Amendola et al. (2018) Amendola, L., Appleby, S., Avgoustidis, A., et al. 2018, Living Reviews in Relativity, 21, 2
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., SipHocz, B. M., et al. 2018, AJ, 156, 123
  * Beck et al. (2020) Beck, R., McPartland, C., Repp, A., Sanders, D., & Szapudi, I. 2020, MNRAS, 493, 2318
  * Bosch et al. (2018) Bosch, J., Armstrong, R., Bickerton, S., et al. 2018, PASJ, 70, S5
  * Doré et al. (2016) Doré, O., Werner, M. W., Ashby, M., et al. 2016, arXiv e-prints, arXiv:1606.07039
  * Doré et al. (2018) Doré, O., Werner, M. W., Ashby, M. L. N., et al. 2018, arXiv e-prints, arXiv:1805.05489
  * Frost et al. (2009) Frost, M. I., Surace, J., Moustakas, L. A., & Krick, J. 2009, ApJ, 698, L68
  * Geach et al. (2017) Geach, J. E., Dunlop, J. S., Halpern, M., et al. 2017, MNRAS, 465, 1789
  * Harikane et al. (2018) Harikane, Y., Ono, Y., Ouchi, M., et al. 2018, PASJ, 70, S11
  * Hartley et al. (2022) Hartley, W. G., Choi, A., Amon, A., et al. 2022, MNRAS, 509, 3547
  * Hildebrandt et al. (2009) Hildebrandt, H., Pielorz, J., Erben, T., et al. 2009, A&A, 498, 725
  * Hyun et al. (2023) Hyun, M., Im, M., Smail, I. R., et al. 2023, ApJS, 264, 19
  * Krick et al. (2008) Krick, J. E., Surace, J. A., Thompson, D., et al. 2008, ApJ, 686, 918
  * Krick et al. (2009) —. 2009, ApJ, 700, 123
  * Laureijs et al. (2011) Laureijs, R., Amiaux, J., Arduini, S., et al. 2011, arXiv e-prints, arXiv:1110.3193
  * Lee et al. (2009) Lee, H. M., Kim, S. J., Im, M., et al. 2009, PASJ, 61, 375
  * Magnier et al. (2016) Magnier, E. A., Schlafly, E. F., Finkbeiner, D. P., et al. 2016, arXiv e-prints, arXiv:1612.05242
  * Magnier et al. (2020a) Magnier, E. A., Sweeney, W. E., Chambers, K. C., et al. 2020a, ApJS, 251, 5
  * Magnier et al. (2020b) Magnier, E. A., Schlafly, E. F., Finkbeiner, D. P., et al. 2020b, ApJS, 251, 6
  * Matsuhara et al. (2006) Matsuhara, H., Wada, T., Matsuura, S., et al. 2006, PASJ, 58, 673
  * Merloni et al. (2012) Merloni, A., Predehl, P., Becker, W., et al. 2012, arXiv e-prints, arXiv:1209.3114
  * Miyazaki et al. (2018) Miyazaki, S., Komiyama, Y., Kawanomoto, S., et al. 2018, PASJ, 70, S1
  * Murakami et al. (2007) Murakami, H., Baba, H., Barthel, P., et al. 2007, PASJ, 59, S369
  * Oi et al. (2021) Oi, N., Goto, T., Matsuhara, H., et al. 2021, MNRAS, 500, 5024
  * Ono et al. (2018) Ono, Y., Ouchi, M., Harikane, Y., et al. 2018, PASJ, 70, S10
  * Radzom et al. (2022) Radzom, B. T., Taylor, A. J., Barger, A. J., & Cowie, L. L. 2022, ApJ, 940, 114
  * Rosenwasser et al. (2022) Rosenwasser, B. E., Taylor, A. J., Barger, A. J., et al. 2022, ApJ, 928, 78
  * Sevilla-Noarbe et al. (2021) Sevilla-Noarbe, I., Bechtol, K., Carrasco Kind, M., et al. 2021, ApJS, 254, 24
  * Shim et al. (2020) Shim, H., Kim, Y., Lee, D., et al. 2020, MNRAS, 498, 5065
  * Sobral et al. (2013) Sobral, D., Smail, I., Best, P. N., et al. 2013, MNRAS, 428, 1128
  * Songaila et al. (2022) Songaila, A., Barger, A. J., Cowie, L. L., Hu, E. M., & Taylor, A. J. 2022, ApJ, 935, 52
  * Songaila et al. (2018) Songaila, A., Hu, E. M., Barger, A. J., et al. 2018, ApJ, 859, 91
  * Taylor et al. (2020) Taylor, A. J., Barger, A. J., Cowie, L. L., Hu, E. M., & Songaila, A. 2020, ApJ, 895, 132
  * Taylor et al. (2021) Taylor, A. J., Cowie, L. L., Barger, A. J., Hu, E. M., & Songaila, A. 2021, ApJ, 914, 79
  * Windhorst et al. (2017) Windhorst, R. A., Alpaslan, M., Ashcraft, T., et al. 2017, JWST Medium-Deep Fields - Windhorst IDS GTO Program, JWST Proposal. Cycle 1, ID. #1176, ,
  * Windhorst et al. (2022) Windhorst, R. A., Cohen, S. H., Jansen, R. A., et al. 2022, arXiv e-prints, arXiv:2209.04119

Here we show the combined catalog columns and descriptions in Table 2.

Table 2: Catalog Columns Column | Unit | Type | Small Catalog | Description
---|---|---|---|---
id |  $\cdots$ | int64 | yes |  hscPipe assigned unique ID number
ra | deg | float64 | yes | Right Ascension (J2000)
dec | deg | float64 | yes | Declination (J2000)
patch |  $\cdots$ | string | yes | Patch id
x | pixel | float32 | yes | x pixel location in mosaic
y | pixel | float32 | yes | y pixel location in mosaic
xx | pixels2 | float32 | yes | xx pixel second moment
yy | pixels2 | float32 | yes | yy pixel second moment
xy | pixels2 | float32 | yes | xy pixel second moment
is_primary |  $\cdots$ | bool | no | Object is not a blended object
num_children |  $\cdots$ | float32 | no | Number of deblended child objects
{$g|r|i|z|y|n816|n921$}_ap2(err) | mag | float32 | yes | {$g|r|i|z|y|n816|n921$} 2$\farcs$diameter aperture magnitude (error)
{$g|r|i|z|y|n816|n921$}_ap3(err) | mag | float32 | no | {$g|r|i|z|y|n816|n921$} 3$\farcs$diameter aperture magnitude (error)
{$g|r|i|z|y|n816|n921$}_ap4(err) | mag | float32 | no | {$g|r|i|z|y|n816|n921$} 4$\farcs$diameter aperture magnitude (error)
{$g|r|i|z|y|n816|n921$}_bkg(err) | mag | float32 | no | {$g|r|i|z|y|n816|n921$} background magnitude (error)
{$g|r|i|z|y|n816|n921$}_kron(err) | mag | float32 | yes | {$g|r|i|z|y|n816|n921$} Kron magnitude (error)
{$g|r|i|z|y|n816|n921$}_cmodel(err) | mag | float32 | no | {$g|r|i|z|y|n816|n921$} CModel magnitude (error)
{$g|r|i|z|y|n816|n921$}_blendedness |  $\cdots$ | float32 | no | {$g|r|i|z|y|n816|n921$} hscPipe Blendedness factor
{$g|r|i|z|y|n816|n921$}_detect |  $\cdots$ | bool | no | Detected in {$g|r|i|z|y|n816|n921$}
{$g|r|i|z|y|n816|n921$}_base_SdssCentroid_flag |  $\cdots$ | bool | no | {$g|r|i|z|y|n816|n921$} centroid fitting error flag
{$g|r|i|z|y|n816|n921$}_base_SdssShape_flag |  $\cdots$ | bool | no | {$g|r|i|z|y|n816|n921$} shape fitting error flag
{$g|r|i|z|y|n816|n921$}_base_PixelFlags_flag_edge |  $\cdots$ | bool | no | Object at edge of {$g|r|i|z|y|n816|n921$} imaging
{$g|r|i|z|y|n816|n921$}_base_PixelFlags_flag_bad |  $\cdots$ | bool | no | Bad pixel in {$g|r|i|z|y|n816|n921$}
{$g|r|i|z|y|n816|n921$}_base_PixelFlags_flag_interpolatedCenter |  $\cdots$ | bool | no | Interpolated pixel in {$g|r|i|z|y|n816|n921$}
{$g|r|i|z|y|n816|n921$}_base_PixelFlags_flag_saturatedCenter |  $\cdots$ | bool | no | Saturated pixel in {$g|r|i|z|y|n816|n921$}
{$g|r|i|z|y|n816|n921$}_base_PixelFlags_flag_crCenter |  $\cdots$ | bool | no | Cosmic ray pixel in {$g|r|i|z|y|n816|n921$}
{$g|r|i|z|y|n816|n921$}_base_PixelFlags_flag_bright_object |  $\cdots$ | bool | no | Near bright object in {$g|r|i|z|y|n816|n921$}
{$g|r|i|z|y|n816|n921$}_modelfit_CModel_flag |  $\cdots$ | bool | no | CModel fitting failed in {$g|r|i|z|y|n816|n921$}


# Rare Charm Decays

M. Destefanis  University of Torino and INFN Torino, Italy, 10125

###### Abstract

Although the Standard Model has been firmly established, the search for
physics beyond the SM is ongoing by investigating new experimental probes.
Rare charm decays are a unique tool to access New Physics studies. The high
luminosity achieved by the modern experiments and their high precision allow
for rare charm decay studies in different scenarios. In this work, some of the
most recent experimental results will be discussed in detail.

##  I Introduction

Three of the four fundamental forces (electromagnetic, weak and strong
interactions) existing in the universe are described by the Standard Model
(SM) theory, which offers as well a classification of all the known particles.
The discovery of the Higgs boson in 2012 allowed to firmly establish the SM.
Despite this, some of the observed phenomena such as baryon-antibaryon
asymmetry, possible existance of dark matter particles, neutrino oscillations
and their non-zero masses, are still unexplained. For this reason, theoretical
and experimental efforts are pointing toward extensions of the SM,
collectively referred to as New Physics (NP) models, in order to define a more
fundamental theory bib:theo1 .

NP searches involve different kind of processes: allowed in SM, forbidden in
SM at tree level, and forbidden in SM bib:bes3white . In the first case, some
processes are expected to hold in SM, but not necessarily in models beyond SM.
In the second case, processes (such as the ones involving flavor-changing
neutral current) which change charm quantum number by one or two units are
forbidden in SM at tree level, but can happen in SM at loop levels, thus
becoming rare. The latter case, includes processes which are allowed by space-
time symmetries, but forbidden in SM, such as baryon or lepton number
violating transitions: they require high statistics, but an observation would
be a hint of physics beyond the SM. The order of the branching fractions of
the SM forbidden processes range from 10-9 to lower values.

##  II Rare Charm Decays

In this scenario, the search for hints of NP was focused on the study of rare
processes in kaon and beauty sectors, in order to find possible discrepancies
with the SM predictions. For long time, rare charm decays have been considered
less promising, due to the different dynamics involved (lack of effective
methods to describe its low energy dynamics) and to the possibly different
couplings with respect to $d$-type quarks bib:brundu . Rare charm decays are
sensitive to $|\Delta c|=|\Delta u|=1$ transitions, but the low statistics
available allows only to set upper limits to the different decay processes. In
order to search for NP, clean null-test observable can be constructed
exploiting exact or approximate simmetries of the SM. Those studies are
complementary to the ones already performed. One of the main issues in this
type of decays is the separation of the short-distance (SD) and long-distance
(LD) information. Indeed, the Glashow–Iliopoulos–Maiani (GIM) mechanism
strongly suppress the Feynman diagrams which describe those decays, so the LD
contributions become dominant. Away from LD contributions, NP should become
visible.

Rare charm decays can be investigated by the modern experiments, due to the
high luminosity and precision achieved. The different experimental scenarios
allow those studies in $e^{+}-e^{-}$ (BESIII bib:bes3 , CLEO-c bib:cleoc ,
BELLE bib:belle , BELLE2 bib:belle2 , and BABAR bib:babar ) and $p-p$ (LHCb
bib:lhcb , CDF bib:cdf , and D0 bib:D0 ) collisions. The first scenario is
carachterized by extremely clean environments, in some cases the double-tag
technique allows almost background free studies, quantum coherence, high
trigger efficiency, and easy detection of neutral particles. On the other
side, the latter scenario offers large production rate, and excellent time
resolution.

In the following, some of the latest results will be discussed.

##  III Flavor Changing Neutral Current

Strong or electromagnetic interactions dominate the $\psi(nS)$ ($n$ = 1, 2)
decays below the open charm threshold. Flavor changing weak decays are anyway
allowed in the SM through the exchange of a virtual $W$ boson. Flavor changing
neutral current (FCNC) transitions $c\rightarrow u\gamma$ and $c\rightarrow
ull$ (where $l$ indicate a lepton) mediate some of the $D$ mesons decays. In
this scenario, the GIM cancellation mechanism strongly suppress those decays;
the suppression in charm decays is much more effective than in kaon or beauty
systems decays due to a stronger diagram cancellation. The SD contributions of
FCNC in the SM are expected to be well beyond the current experimental
sensitivity, but some theoretical model suggests that they could be enhanced
by LD effects by several orders of magnitude. Hence, larger FCNC transition
rates are predicted by some NP scenario: the observation of such an
enhancement would directly indicate NP.

###  III.1 $D^{0}\rightarrow h^{+}h^{-}\mu^{+}\mu^{-}$

The LHCb Collaboration, taking advantage of p-p annihilation data collected
during Run 1 (7 and 8 TeV) and Run 2 (13 TeV) for a total of 9fb-1,
investigated the decay $D^{*+}\rightarrow D^{0}\pi^{+}$, $D^{0}\rightarrow
h^{+}h^{-}\mu^{+}\mu^{-}$, where $h$ indicate either a pion or a kaon
bib:D0hhmumu . The semileptonic decays are described by means of five
independent kinematic variables: $q^{2}=m^{2}(\mu^{+}\mu^{-})$ and
$p^{2}=m^{2}(h^{+}h-)$, the dimuon and dihadron squared invariant masses, and
the three decay angles $\vartheta_{\mu}$, $\vartheta_{h}$, and $\varphi$,
defined as depicted in Fig. 1.

Figure 1: Decay topology of $D^{0}\rightarrow h^{+}h^{-}\mu^{+}\mu^{-}$.

The differential decay rate can be described as the sum of nine angular
coefficients $I_{1-9}$, which depend on $q^{2}$, $p^{2}$, and
$\cos{\vartheta_{h}}$, multiplied by $c_{1-9}$ terms, functions of
$\vartheta_{\mu}$ and $\varphi$. $I_{1}$ is a normalization factor, while
$I_{2-9}$ are expressed as angular asymmetries. The decay rate asymmetries are
the observables experimentally measured, for example
$<I_{2}>=\frac{\Gamma(|cos\vartheta_{\mu}|>0.5)-\Gamma(|cos\vartheta_{\mu}|<0.5)}{\Gamma(|cos\vartheta_{\mu}|>0.5)+\Gamma(|cos\vartheta_{\mu}|<0.5)}$,
where $\Gamma$ is the decay rate. For those angular asymmetries it is possible
to define their $CP$ average as
$<S_{i}>=\frac{1}{2}[<I_{i}>+(-)<\bar{I_{i}}>]$ and the asymmetry as
$<A_{i}>=\frac{1}{2}[<I_{i}>-(+)<\bar{I_{i}}>]$. The observables $<S_{5,6,7}>$
and the $CP$ asymmetries $<A_{2-9}>$ are expected to vanish and they can
constitute the SM null test, if only SM amplitudes contribute to the decay
processes. The D0 invariant masses were fitted with an Hypatia distribution
for the signal and a Johnson $S_{U}$ distribution plus an exponential function
for the background. The signal yields for
$D^{0}\rightarrow\pi^{+}\pi^{-}\mu^{+}\mu^{-}$ and $D^{0}\rightarrow
K^{+}K^{-}\mu^{+}\mu^{-}$ were found to be 3579$\pm$71 and 318$\pm$19,
respectively. The values for the observables $<S_{5,6,7}>$ and $<A_{2-9}>$ are
in agreement with the SM null hypothesis.

###  III.2 $D^{0}\rightarrow\pi^{0}\nu\bar{\nu}$

The BESIII Collaboration searched for neutrino pair in the final state of
charm decays for the first time bib:PRD105_L071102 , taking advantage of 2.93
fb-1 of data collected at the energy of 3.773 GeV. The $e^{+}e^{-}\rightarrow
D^{0}\bar{D}^{0}$ process, where $D_{0}\rightarrow\pi^{0}\nu\bar{\nu}$ was
investigated. In this decay, LD contributions become insignificant and the SD
contributions from Z-penguin and box diagrams are dominant. The final state
was reconstructed by means of the double-tag technique: first the single-tag
$\bar{D}^{0}$ was reconstructed in three different hadronic channels
($K^{+}\pi^{-}$, $K^{+}\pi^{-}\pi^{0}$, and $K^{+}\pi^{-}\pi^{+}\pi^{-}$),
then the $D^{0}$ decay signal was reconstructed recoiling against the single-
tag $\bar{D}^{0}$, searching for two $\gamma$s from $\pi^{0}$ decay. Two more
variables were investigated: the beam-constrained mass
$M_{BC}=\sqrt{E^{2}_{beam}/c^{4}-|\bf{p}_{\bar{D}^{0}}|^{2}/c^{2}}$ and the
energy difference $\Delta E=E_{\bar{D}^{0}}-E_{beam}$, where $E_{beam}$ is the
energy of the electron beam. In order to be able to investigate the
$D_{0}\rightarrow\pi^{0}\nu\bar{\nu}$ decay, a reliable modeling of the
background contributions is mandatory. To achieve this result, a data-driven
method was used: the data-driven analysis procedure was validated by means of
one third of the full data sample. Finally, the signal was searched in the
energy deposited in the electromagnetic calorimeter. The energy distribution,
shown in Fig. 2, depict the data distribution with a global fit and the
different contributions. The signal, modeled on MC simulations, is indicated
by the gray area and it is normalized to 20 times the central value of the fit
result for visibility; 14$\pm$30 signal events were found. Since the number of
reconstructed events is compatible with zero, the branching ratio upper limit
for the $D_{0}\rightarrow\pi^{0}\nu\bar{\nu}$ decay was calculated to be
2.1$\times$10-4 at 90% confidence level (C.L.).

Figure 2: Fit of the energy deposited in the electromagnetic calorimeter and
different contributions together with the pull. Dots with error bars indicate
the data, the red solid line indicates the global fit, the green long dashed
line indicates the $D^{0}\rightarrow\pi^{0}K^{0}_{L}X$ contribution, the blue
middle dashed line include the other $D^{0}$ decays, the magenta dashed line
indicates the wrong tags, and the signal is depicted by the grey shaded area.

##  IV Lepton Flavor Violating processes

In the charm sector, lepton flavor conserving processes, such as $D\rightarrow
Xe^{+}e^{-}$ and $D\rightarrow X\mu^{+}\mu^{-}$, are predicted in the SM. They
can occurr through SD and LD processes and they are expected to have branching
fractions of the order of $\mathcal{O}(10^{-9})$ and $\mathcal{O}(10^{-6})$,
respectively bib:PRD101_112003 . The lepton flavor violating (LFV) neutral
decays, such as $D^{0}\rightarrow X^{0}e^{\pm}\mu^{\mp}$, are prohibited in
the SM, since they can occurr only through lepton mixing bib:PLB751_54 ; their
branching ratio are expected to be in the order of $\mathcal{O}(10^{-50})$.
However, different NP models, such as leptoquarks, two-Higgs doublets, and
those involving Majorana neutrinos, allow both LFV and lepton number violation
bib:LFVpapers , with branching fractions of the order of
$\mathcal{O}(10^{-5})$. The discovery of neutrino oscillations confirmed the
LFV in the neutral lepton sector as well as the existence of neutrino masses,
while a detection of LFV in the charged lepton sector would provide direct
evidence of NP bib:bes3white .

###  IV.1 $D^{0}\rightarrow X^{0}e^{\pm}\mu^{\mp}$

The BABAR Collaboration investigated the decay $D^{*+}\rightarrow
D^{0}\pi^{+}$, $D^{0}\rightarrow X^{0}e^{\pm}\mu^{\mp}$, where $X^{0}$
indicates a neutral meson ($\pi^{0}$, $K_{S}^{0}$, $\bar{K}^{*0}$, $\rho^{0}$,
$\varphi$, $\omega$, or $\eta$), exploiting 424 fb-1 of data collected at the
$\Upsilon(4S)$ resonance (10.58 GeV) and additional 44 fb-1 of data collected
0.04 GeV below the $\Upsilon(4S)$ resonance bib:PRD101_112003 . The decays
$D^{0}\rightarrow\pi^{-}\pi^{+}\pi^{+}\pi^{-}$ (for
$X^{0}=K_{S}^{0},\rho^{0},\omega$), $D^{0}\rightarrow
K^{-}\pi^{+}\pi^{+}\pi^{-}$ (for $X^{0}=\bar{K}^{*0}$), and $D^{0}\rightarrow
K^{-}K^{+}\pi^{+}\pi^{-}$ (for $X^{0}=\varphi$) were used to normalize the
branching fractions for the different signal modes. After the full
reconstruction of the signal modes, an unbinned maximum-likelyhood fit is
performed to the variable $\Delta m=m(D^{*+})-m(D^{0})$. The distributions for
each signal mode are shown in Fig. 3, where dots with error bars are the
experimental data, the red dashed line is the signal, the green dotted line is
the background, end the black solid line is the total fit. No significant
signal was observed; the extracted upper limits at 90% C.L. are reported in
Tab. 1. The obtained values are between 1 and 2 orders of magnitude more
stringent with respect to the previous results.

Figure 3: Unbinned maximum-likelyhood fits to the final candidate distributions as function of $\Delta m$. Table 1: Summary of branching ratios upper limits at 90% C.L. for the investigated decay modes. Decay mode | BR U.L. ($\times 10^{-7}$)
---|---
$D^{0}\rightarrow\pi^{0}e^{\pm}\mu^{\mp}$ | 8.0
$D^{0}\rightarrow K^{0}_{S}e^{\pm}\mu^{\mp}$ | 8.7
$D^{0}\rightarrow\bar{K}^{*0}e^{\pm}\mu^{\mp}$ | 12.5
$D^{0}\rightarrow\rho^{0}e^{\pm}\mu^{\mp}$ | 5.0
$D^{0}\rightarrow\varphi e^{\pm}\mu^{\mp}$ | 5.1
$D^{0}\rightarrow\omega e^{\pm}\mu^{\mp}$ | 17.1
$D^{0}\rightarrow\eta e^{\pm}\mu^{\mp}$ | 22.5
with $\eta\rightarrow\gamma\gamma$ | 24.0
with $\eta\rightarrow\pi^{+}\pi^{-}\pi^{0}$ | 43.0

##  V Lepton Number Violating processes

In the SM, neutrinos are postulated to be massless, due to the absence of
right-handed neutrino. However, the first evidence for physics beyond the SM
is due to the small neutrino masses observed as consequence of neutrino
oscillation. Moreover, it is still an open question whether neutrinos are
Dirac or Majorana particles. If they are Majorana particles, lepton number
violating (LNV) processes by two units ($\Delta L=2$) can be observed.
Different NP models involving LNV have been proposed at different energy
regime bib:LNVpapers . The exchange of a single Majorana neutrino with a mass
on the order of the heavy flavor mass scale could be a source of LNV
processes.

###  V.1 $D\rightarrow K\pi e^{+}e^{+}$

The BESIII Collaboration searched for LNV processes with $\Delta L=2$ in D
meson decays bib:PRD99_112002 , taking advantage of 2.93 fb-1 of data
collected at the energy of 3.773 GeV. The processes $D^{0}\rightarrow
K^{-}\pi^{-}e^{+}e^{+}$, $D^{+}\rightarrow K^{0}_{S}\pi^{-}e^{+}e^{+}$, and
$D^{+}\rightarrow K^{-}\pi^{0}e^{+}e^{+}$ are mediated by a Majorana neutrino
$\nu_{m}$ and can occur through Cabibbo-favoured (CF) and doubly Cabibbo-
suppressed (DCS) decays. The latter are expected to be suppressed by a factor
0.05 with respect to CF processes. In this analysis, the single-tag method was
used. For the signal determination, two variables were taken into account: the
beam energy constrained mass $M_{BC}=\sqrt{E^{2}_{beam}-|\bar{p}_{D}|^{2}}$
and the energy difference $\Delta E=E_{D}-E_{beam}$, where $E_{beam}$ is the
beam energy, and $E_{D}$ and $\bar{p}_{D}$ are the energy and momentum of $D$
candidates. In order to extract the signal yields, an unbinned maximum
likelyhood fit was performed on the $M_{BC}$ distributions; the signal was
described by a MC simulated shape convolved with a Gaussian, while the
background was described by an ARGUS function. No obvious signal was observed.
By means of the Bayesian approach, the branching fraction upper limits at 90%
CL were calculated to be 2.8$\times$10-6 for $D^{0}\rightarrow
K^{-}\pi^{-}e^{+}e^{+}$, 3.3$\times$10-6 for $D^{+}\rightarrow
K^{0}_{S}\pi^{-}e^{+}e^{+}$, and 8.5$\times$10-6 for $D^{+}\rightarrow
K^{-}\pi^{0}e^{+}e^{+}$. With different $\nu_{m}$ mass hypothesis the
branching fractions upper limits at 90% CL for the processes $D^{0}\rightarrow
K^{-}e^{+}\nu_{m}(\pi^{-}e^{+})$ and $D^{+}\rightarrow
K^{0}_{S}e^{+}\nu_{m}(\pi^{-}e^{+})$ were found to be at level of 10-6-10-7,
as shown in Fig. 4.

Figure 4: Upper limits on branching fractions as function of different
$\nu_{m}$ masses for the decays (a) $D^{0}\rightarrow
K^{-}e^{+}\nu_{m}(\pi^{-}e^{+})$ and (b) $D^{+}\rightarrow
K^{0}_{S}e^{+}\nu_{m}(\pi^{-}e^{+})$.

##  VI Baryon Number Violating processes

In the SM, as a consequence of the $SU(2)\times U(1)$ and $SU(3)$ gauge
symmetries, the baryon number is conserved. However, the existence of baryon
number violating (BNV) processes is suggested by the excess on baryon with
respect to antibaryon in the Universe. Thus, the investigation of BNV
processes could allow to understand the evoulution of the Universe. Some
extensions of the SM include BNV processes bib:BNVpapers , for example under
dimension six operators BNV can happen with a change in the baryon and lepton
numbers $\Delta(B-L)=0$, or under dimension seven operators $\Delta(B-L)=2$ is
allowed. The Feynman diagrams describing those processes include the presence
of heavy gauge boson $X$ with charge $\frac{4}{3}$ and gauge boson $Y$ with
charge $\frac{1}{3}$, which can couple a quark to a lepton.

###  VI.1 $D^{+}\rightarrow\bar{\Lambda}(\bar{\Sigma}^{0})e^{+}$

The BESIII Collaboration measured for the first time BNV processes with
$\Delta(B-L)=0$ in $D^{+}\rightarrow\bar{\Lambda}(\bar{\Sigma}^{0})e^{+}$ and
$\Delta(B-L)=2$ in $D^{+}\rightarrow\Lambda(\Sigma^{0})e^{+}$ decays
bib:PRD101_031102 , taking advantage of 2.93 fb-1 of data collected at the
energy of 3.773 GeV. The two hyperons were reconstructed through the channels
$\Lambda\rightarrow p\pi^{-}$ and $\Sigma^{0}\rightarrow\gamma\Lambda$. The
BNV decays were investigated by means of the beam energy constrained mass
$M_{BC}=\sqrt{E^{2}_{beam}-|\bar{p}_{D}|^{2}}$ and the energy difference
$\Delta E=E_{D}-E_{beam}$, where $E_{beam}$ is the beam energy, and $E_{D}$
and $\bar{p}_{D}$ are the energy and momentum of $D^{+}$ candidates. Figure 5
shows the performed maximum likelyhhod fit to the $M_{BC}$ distributions; here
the signal is modeled with a MC simulated shape convolved with a Gaussian,
while the background is described by an ARGUS function. No obvious signal was
found. The branching fractions upper limits at 90% CL were found to be
$1.1\times 10^{-6}$ for $\Lambda e^{+}$, $6.5\times 10^{-7}$ for
$\bar{\Lambda}e^{+}$, $1.7\times 10^{-6}$ for $\Sigma^{0}e^{+}$, and
$1.3\times 10^{-6}$ for $\bar{\Sigma}^{0}e^{+}$.

Figure 5: Fits to the $M_{BC}$ distributions. Dots with error bars are the
data, the solid lines are the fits, the red dashed lines indicate the
background, and the blue hatched histograms are the MC simulated backgrounds
scaled to the data according to the luminosity bib:PRD101_031102 .

###  VI.2 $D^{0}\rightarrow pe^{-}$

The BESIII Collaboration measured baryon and lepton production in $D^{0}$
meson decays bib:PRD105_032006 , taking advantage of 2.93 fb-1 of data
collected at the energy of 3.773 GeV. The double-tag technique was implied to
reconstruct the $\psi(3770)$ decays into $D^{0}\bar{D}^{0}$ meson pairs. At
first, the $\bar{D}^{0}$ was reconstructed via its hadronic decay modes
$K^{+}\pi^{-}$, $K^{+}\pi^{-}\pi^{0}$, and $K^{+}\pi^{-}\pi^{-}\pi^{+}$
(single-tag ST). Two more variables were investigated: the beam-constrained
mass $M_{BC}=\sqrt{E^{2}_{beam}/c^{4}-|\bf{p}_{\bar{D}^{0}}|^{2}/c^{2}}$ and
the energy difference $\Delta E=E_{\bar{D}^{0}}-E_{beam}$, where $E_{beam}$ is
the energy of the electron beam, and $E_{\bar{D}^{0}}$ and ${\bf
p}_{\bar{D}^{0}}$ are the energy and the momentum of the candidate
$\bar{D}^{0}$. In order to extract the yield of the ST mesons, the $M_{BC}$
distributions were fitted with a MC simulated shape convolved with a double-
Gaussian for the signal and an ARGUS function for the backgroud;
2321009$\pm$1875 events were reconstructed. In order to extract the double-tag
signal $D^{0}\rightarrow\bar{p}e^{+}$ and $D^{0}\rightarrow pe^{-}$, the
variables $M_{BC}^{sig}$ and $\Delta E^{sig}$ were defined in a similar way to
the ST ones. The signal yields were obtained by counting the events in the
selected region as shown in Fig. 6; no obvious signals were found. The
branching fractions upper limits at 90% CL were found to be $<1.2\times
10^{-6}$ for $D^{0}\rightarrow\bar{p}e^{+}$, and $<2.2\times 10^{-6}$ for
$D^{0}\rightarrow pe^{-}$.

Figure 6: $M_{BC}^{sig}$ as function of $\Delta E^{sig}$ for $D^{0}\rightarrow
pe^{-}$. The red rectangle indicates the signal region.

###  VI.3 $D_{(S)}^{+}\rightarrow hll$

The LHCb Collaboration, taking advantage of p-p annihilation data collected
during Run 1 at 8 TeV for a total of 1.6fb-1, investigated the decays
$D^{+}_{(S)}\rightarrow h^{\pm}l^{+}l^{(^{\prime})\mp}$, where $h$ indicate
either a pion or a kaon, and $l^{(^{\prime})\mp}$ is an electron or a muon
bib:jhep06_044 . FCNC transitions are involved in four decay channels
($D^{+}\rightarrow\pi^{+}e^{+}e^{-},\pi^{+}\mu^{+}\mu^{-}$, and
$D^{+}_{S}\rightarrow K^{+}e^{+}e^{-},K^{+}\mu^{+}\mu^{-}$); LFV processes are
involved in eight decays
($D^{+}_{(S)}\rightarrow\pi^{+}e^{+}\mu^{-},\pi^{+}\mu^{+}e^{-},K^{+}e^{+}\mu^{-},K^{+}\mu^{+}e^{-}$).
Nine decay modes involve both LFV and LNV
($D^{+}_{(S)}\rightarrow\pi^{-}e^{+}e^{+},\pi^{-}\mu^{+}\mu^{+},\pi^{-}\mu^{+}e^{+}$,
and $D^{+}_{S}\rightarrow
K^{-}e^{+}e^{+},K^{-}\mu^{+}\mu^{+},K^{-}\mu^{+}e^{+}$). For calibration and
normalization four resonant decays (NRD)
($D^{+}_{(S)}\rightarrow(\varphi\rightarrow\mu^{-}\mu^{+})\pi^{+},(\varphi\rightarrow
e^{-}e^{+})\pi^{+}$) were exploited. Those decays are dominated by LD tree
level contributions, and the regions dominated by resonances in dilepton mass
were vetoed when fitting. The NRD invariant masses were fitted with a double
Gaussian and an exponential function or a third order Chebyshev polynomial for
muons and electrons final states, respectively. No significant deviation from
the background was found for the 25 decay modes investigated. The branching
fractions upper limits at 90% CL are reported in Tab. 2, and represent an
improvement of 1 or 2 orders of magnitude with respect to the existing limits.

Table 2: Upper limits at 90% CL for each signal decay channel $D_{(S)}^{0}\rightarrow hll$. Decay mode | BR U.L. ($\times 10^{-9}$)
---|---
$D^{+}$ |  $D^{+}_{S}$
$D^{0}_{(S)}\rightarrow\pi^{+}\mu^{+}\mu^{-}$ |  67 |  180
$D^{0}_{(S)}\rightarrow\pi^{-}\mu^{+}\mu+$ |  14 |  86
$D^{0}_{(S)}\rightarrow K^{+}\mu^{+}\mu^{-}$ |  54 |  140
$D^{0}_{(S)}\rightarrow K^{-}\mu^{+}\mu+$ |  - |  26
$D^{0}_{(S)}\rightarrow\pi^{+}e^{+}\mu^{-}$ |  210 |  1100
$D^{0}_{(S)}\rightarrow\pi^{+}\mu^{+}e^{-}$ |  220 |  940
$D^{0}_{(S)}\rightarrow\pi^{-}\mu^{+}e^{+}$ |  130 |  630
$D^{0}_{(S)}\rightarrow K^{+}e^{+}\mu^{-}$ |  75 |  790
$D^{0}_{(S)}\rightarrow K^{+}\mu^{+}e^{-}$ |  100 |  560
$D^{0}_{(S)}\rightarrow K^{-}\mu^{+}e^{+}$ |  - |  260
$D^{0}_{(S)}\rightarrow\pi^{+}e^{+}e^{-}$ |  1600 |  5500
$D^{0}_{(S)}\rightarrow\pi^{-}e^{+}e^{+}$ |  530 |  1400
$D^{0}_{(S)}\rightarrow K^{+}e^{+}e^{-}$ |  850 |  4900
$D^{0}_{(S)}\rightarrow K^{-}e^{+}e^{+}$ |  - |  770

###  VI.4 $D^{0}\rightarrow hhll$

The BABAR Collaboration investigated the decay $D^{*+}\rightarrow
D^{0}\pi^{+}$, with $D^{0}$ decays into $h^{\prime-}h^{-}l^{\prime+}l^{+}$
(nine LNV dacays) or into $h^{\prime-}h^{+}l^{\prime\pm}l^{\mp}$ (three LFV
decays), where $h$ and $h^{\prime}$ indicate a K or a $\pi$ meson, and $l$ and
$l^{\prime}$ an electron or a muon bib:PRL124_071802 . The 424 fb-1 of data
collected at the $\Upsilon(4S)$ resonance (10.58 GeV) and additional 44 fb-1
of data collected 0.04 GeV below the $\Upsilon(4S)$ resonance were exploited.
The decays $D^{0}\rightarrow\pi^{-}\pi^{+}\pi^{+}\pi^{-}$, $D^{0}\rightarrow
K-\pi^{+}\pi^{+}\pi^{-}$, and $D^{0}\rightarrow K^{-}K^{+}\pi^{+}\pi^{-}$ were
used to normalize the branching fractions for the different signal modes. In
order to reject semileptonic charm decays or final states containing a neutral
particle, a multivariate selection, based on a Fisher discriminant which uses
nine input observables, was applied; the discriminants were trained and tested
using MC for the signal modes. An unbinned maximum-likelyhood fit was
performed to the variable $\Delta m=m(D^{*+})-m(D^{0})$ for the normalization
modes and for the signal ones. In the fit, the signal was described with the
sum of multiple Cruijff and Crystal Ball functions, while the background was
described with an ARGUS function. For all the signal modes, the yields were
found compatible with zero. The extracted upper limits at 90% C.L. are
reported in Tab. 3. The obtained values are between 1 and 3 order of magnitude
more stringent with respect to the previous results.

Table 3: Upper limits at 90% CL for each signal decay channel $D^{0}\rightarrow hhll$. $D^{0}$ decay mode | BR U.L. ($\times 10^{-7}$)
---|---
$\pi^{-}\pi^{-}e^{+}e^{+}$ | 9.1
$\pi^{-}\pi^{-}\mu^{+}\mu+$ | 15.2
$\pi^{-}\pi^{-}e^{+}\mu^{+}$ | 30.6
$\pi^{-}\pi^{+}e^{\pm}\mu^{\mp}$ | 17.1
$K^{-}\pi^{-}e^{+}e^{+}$ | 5.0
$K^{-}\pi^{-}\mu^{+}\mu+$ | 5.3
$K^{-}\pi^{-}e^{+}\mu^{+}$ | 21.0
$K^{-}\pi^{+}e^{\pm}\mu^{\mp}$ | 19.0
$K^{-}K^{-}e^{+}e^{+}$ | 3.4
$K^{-}K^{-}\mu^{+}\mu+$ | 1.0
$K^{-}K^{-}e^{+}\mu^{+}$ | 5.8
$K^{-}K^{+}e^{\pm}\mu^{\mp}$ | 10.0

##  VII Future Plans

The BelleII experiment will have the capability to investigate rare or
forbidden charm decays. In particular, those involving neutral particles or
missing energy in the the final states (such as neutrinos, dark matter,
axions, other non-SM particles) will be well suited. Among the others, the
$D^{0}\rightarrow\gamma\gamma$, $c\rightarrow u$ FCNC will be unique for the
experiment. The branching fraction upper limit for this process is about
$8.7\times 10^{-7}$ bib:belle_D0gg  and it is still two order of magnitude
higher with respect to the SM predictions; the expactation is to reach about 1
order of magnitude more stringent upper limit, when 50 ab-1 of data will be
collected. The BelleII Collaboration plans as well to study the $D\rightarrow
h\nu\bar{\nu}$ and the $\Lambda_{c}^{+}\rightarrow p\nu\bar{\nu}$ final
states, which are strongly suppressed in the SM due to the GIM mechanism,
taking advantage of the powerful reconstruction method developed.

On the other hand, the BESIII Collaboration, which get an extended operation
time, is planning to collect more data in the charm region. In particular, 20
fb-1 and 6 fb-1 are expected at the 3.773 GeV and 4.18 GeV energies, in order
to further improve our knowledge on rare charm decays bib:bes3white .

##  VIII Conclusions

Rare charm decays are a great tool for NP investigations. They can pin down
the $u$-quark dynamics, and deliver complementary information with respect to
the strange and the beauty sectors. Unfortunately, the collected statistics is
increasing, but it is still too low to be able to reach the SM predictions. So
far, no NP effect has been found. New results are expected from the different
experiments and, possibly, from next generation facilities.

###### Acknowledgements.

I would like to thank the BESIII, BelleII, LHCb, and BABAR Collaborations for
providing me the excellent results discussed in this paper.

## References

  * (1) H. Gisbert et al., Mod. Phys. Lett. A 36, No. 4, 2130002 (2021).
  * (2) M. Ablikim et al., Chin. Phys. C 44, 040001 (2020).
  * (3) D. Brundu, Il Nuovo Cimento 42 C, 159 (2019).
  * (4) D.M. Asner et al., arXiv:0809.1869v1 [hep-ex] (2008).
  * (5) R.A. Briere et al., CLNS 01/1742, Revised 10/01 (2001).
  * (6) A. Abashian et al., Nucl. Instrum. Meth. A 479, 117 (2002).
  * (7) T. Abe et al., arXiv:1011.0352 (2010); E. Kou et al., Prog. Theor. Exp. Phys. 2019, 123C01 (2019).
  * (8) B. Aubert et al., Nucl. Instrum. Meth. A 479, 1-116 (2002).
  * (9) The LHCb Collaboration et al., JINST 3, S08005 (2008).
  * (10) F. Abe et al., Nucl. Instrum. Meth. A 271, 387-403 (1998).
  * (11) S. Abachi et al., Nucl. Instrum. Meth. A 338, 185-253 (1994).
  * (12) R. Aaij et al., Phys. Rev. Lett. 128, 221801 (2022).
  * (13) M. Ablikim et al., Phys. Rev. D 105, L071102 (2022).
  * (14) D. Guadagnoli and K. Lane, Phys. Lett. B 751, 54 (2015).
  * (15) S. de Boer and G. Hiller, Phys. Rev. D 93, 074001 (2016); S. Fajfer and N. Kosnik, Eur. Phys. J. C 75, 567 (2015); A. Atre et al., J. High Energy Phys. 05, 030 (2009); H. Yuan et al., J. High Energy Phys. 08, 066 (2013); H.R. Dong et al., Chin. Phys. C 39, 013101 (2015).
  * (16) J.P. Lees et al., Phys. Rev. D 101, 112003 (2020).
  * (17) A. Atre et al., J. High Energy Phys. 05, 030 (2009); G.L. Castro and N. Quintero, Phys. Rev. D 87, 077901 (2013); H.R. Dong et al., Chin. Phys. C 39, 013101 (2015); H.B. Li, Front. Phys. 12, 121301 (2017); D. Milanés and N. Quintero, Phys. Rev. D 98, 096004 (2018).
  * (18) M. Ablikim et al., Phys. Rev. D 99, 112002 (2019).
  * (19) S. Weinberg, Phys. Rev. Lett. 43, 1566 (1979); S. Weinberg, Phys. Rev. D 22, 1694 (1980); F. Wilczek and A. Zee, Phys. Rev. Lett. 43, 1571 (1979); H.A. Weldon and A. Zee, Nucl. Phys. 173, 269 (1980); F. Wilczek and A. Zee, Phys. Lett. B 88, 311 (1979); R.N. Mohapatra and R.E. Marshak, Phys. Lett. B 94, 183 (1980); K. Biswal et al., Phys. Rev. D 25, 266 (1982); W.S. How et al.,Phys. Rev. D 72, 095001 (2005).
  * (20) M. Ablikim et al., Phys. Rev. D 101, 031102 (2020).
  * (21) M. Ablikim et al., Phys. Rev. D 105, 032006 (2022).
  * (22) R. Aaij et al., JHEP06(2021)044.
  * (23) J.P. Lees et al., Phys. Rev. Lett. 124, 071802 (2020).
  * (24) N.K. Nisar et al., Phys. Rev. D 93, 051102 (2016).