things: vortexes and various structures. Philosophers can discuss if things
exist or are created by us in the course of cognition. But the practical point
of view is undoubted: the hypothesis of things is very useful and there is an
enormous amount of information in it about the world. (In some sense, the
atomic hypothesis can be read as an extension of the hypothesis of things to
the microscale.)

We can go even further and apply this hypothesis of things to features of
objects and to relations between them, like L. Boltzmann advised to M. Planck
to apply the atomic hypothesis to emission of electromagnetic radiation (and
we remember the great success that followed this advice). We will find that
features and relations are also combined into distinct discrete entities. The
main purpose of big data preprocessing is to uncover all of these entities
hidden in the data. After that, various pattern recognition tasks become
viable.

The preprocessing in the postclassical data world (Sec. 2, Fig. 3) is a
challenging task because no classical statistical methods are applicable when
the sample size is much less than data dimensionality (the Donoho area (Sec.
2.1, 1) Donoho (2000)). The correlation transformation (Sec. 2.2) moves data
out of the Donoho area, but the specific non-classical effects persist while
the sample size remains much less than the exponential of the data
dimensionality (2). Dimensionality reduction methods should combine two set of
goals: sensible grouping and extraction of relevant features. For these
purposes, combination of supervised and unsupervised learning techniques is
necessary. Data labels from supervised approaches add sense to analysis of
unlabeled data. The simple geometric methods like supervised PCA,
semisupervised PCA (Sec. 2.3), and Domain Adaptation PCA (DAPCA) (Sec. 2.3)
serve as prototypes of more complex and less controllable approaches. They can
also be used to simplify large deep learning systems GorbanMirkesTukin2020 .

Data in postlclassical world are rarefied. At the same time, values of regular
functionals on data are concentrated near their median values Giannopoulos and
Milman (2000); Vershynin (2018). Combinations of these properties produce
‘blessing of dimensionality’ Donoho (2000); Kainen1997 (1997); Lévy (1951).
The most important manifestation of these effects for applied data analysis
beyond the central limit theorem are quasiorthogonality Kainen & Kůrková
(1993, 2020); Gorban et al. (2016a) and stochastic separation theorems
GorbanTyukinNN2017 ; Grechuk2021 . These results give the theoretical
backgrounds for construction of new type of intellectual devices, correctors
of AI systems. In this paper, we presented a new family of stochastic
separation theorems for fine-grained data distributions with different
geometry of clusters (Sec. 3). These results allowed us to develop the multi-
correctors of multidimensional AI with a granular distribution of errors. On
real data, such correctors showed better performance than simple correctors.

Various versions of multi-correctors that provide fast and reversible
correction of AI errors should be supplemented by the special operation of
interiorization of corrections. Accumulation of many corrections will. step by
step, spend the blessing of dimensionality resource: after implementing
elementary corrections, the probability of success of new correctors
decreases, which can be considered as accumulation of technical debt. In
psychology, interiorization is the process of making skills, attitudes,
thoughts, and knowledge an integrated part of one’s own being. For large
legacy AI systems, interiorization of corrections means the supervising re-
training of the system, where the complex “legacy system+multi-corrector” acts
as a supervisor and labels the data, while the system itself learns by
assimilating the fast flow of generated data.

The infinite-dimensional version of theorems about separation of compact
clusters and families of such clusters demonstrates the importance of
hypothesis about compact embedding of data clusters (Sec. 3.5). The idealized
concept of granular Hierarchical Universe (Fig. 5) is intended to replace the
ideal picture of a smooth unimodal distribution popular in statistical
science.

In the future, the concept of granular Hierarchical Universe will be extended
to feature space and to relationships between different entities (in data
space, feature space, etc.), which is natural in the deep learning framework.

conceptualization and methodology, ANG and IYT; writing–original draft
preparation, ANG, IYT and BG; writing–review and editing, all authors;
software and validation, IYT and EMM.

IYT was funded by UKRI (Alan Turing AI Acceleration Fellowship EP/V025295/1).
ANG, EMM, and IYT were founded by the Ministry of Science and Higher Education
of the Russian Federation (Project No. 075-15-2020-808).

The authors declare no conflict of interest. The funders had no role in the
design of the study; in the collection, analyses, or interpretation of data;
in the writing of the manuscript, or in the decision to publish the results.

The following abbreviations are used in this manuscript:

AI | Artificial Intelligence
---|---
ML | Machine Learning
PCA | Principal Component Analysis
TCA | Transfer Component Analysis
DAPCA | Domain Adaptation PCA

##  Appendix A  “Almost Always" in infinite-dimensional spaces

As it was mentioned in Sec. 3.5, in the infinite-dimensional limit many
statements about high or low probabilities transform into 0-1 laws: something
happens almost always or almost newer. Such limits for concentrations on
spheres and their equators were discussed by Lévy Lévy (1951) as an important
part of the measure concentration effects. In physics, this limit corresponds
to the so-called thermodynamic limit of statistical mechanics Khinchin1949 ;
Thompson2015 . The original Kolmogorov 0-1 law states, roughly speaking, that
an event that depends on an infinite collection of independent random
variables but is independent of any finite subset of these variables, has
probability zero or one (for precise formulation we refer to the monograph
Kolmogorov2018 ). The infinite-dimensional 0-1 asymptotic might bring more
light and be more transparent than the probabilistic formulas.

This may be surprising, but the problem is what “almost always” means.
Formally, various definitions of genericity are constructed as follows. All
systems (or cases, or situations and so on) under consideration are somehow
parameterized – by sets of vectors, functions, matrices etc. Thus, the “space
of systems" $Q$ can be described. Then the “meagre (or thin) sets" are
introduced into $Q$, i.e. the sets, which we shall later neglect. The union of
a finite or countable number of meager sets, as well as the intersection of
any number of them should be meager set again, while the whole $Q$ is not
thin. There are two traditional ways to determine thinness.

  1. 1.

The sets of measure zero are negligible.

  2. 2.

The sets of Baire first category are negligible.

The first definition requires existence of a special measure such that all
relevant distributions are expected to be absolute continuous with respect to
it. In Theorem 3.1, for example, we assumed that the probability distribution
(yet unknown) has density and is absolutely continuous with respect Lebesgue
measure. Moreover, we used a version of the “Smeared (or Smoothed) Absolute
Continuity” (SmAC) condition (11) GorbanGrechukTykin2018 ; Gorban et al.
(2018), which means that the sets of relatively small volume cannot have high
probability, whereas absolute continuity means that sets of zero volume have
probability zero. Unfortunately, in the infinite-dimensional spaces we usually
do not have such a sensible measure. It is very easy to understand if we look
on the volumes of balls in Hilbert space with orthonormal basis
$\\{\boldsymbol{e}_{i}\\}$. If the measure of a ball is function of its radius
and the measure of a ball of radius $R$ is finite, then the balls of radius
$R/4$ have zero measure (because infinitely many such balls with the centers
at points $R\boldsymbol{e}_{i}/2$ can be packed in the ball of radius $R/4$),
and, therefore, the ball of radius $R$ has zero measure because it can be
covered by a countable set of balls of radius $R/4$. Hence, all balls have
either zero or infinite measure.

The second definition is widely accepted when we deal with the functional
parameters. The construction begins with nowhere dense sets. The set $Y$ is
nowhere dense in $Q$, if in any nonempty open set $V\subset Q$ (for example,
in a ball) there exists a nonempty open subset $W\subset V$ (for example, a
ball), which does not intersect with $Y$: $W\cap Y=\emptyset$. Roughly
speaking, Y is “full of holes" – in any neighborhood of any point of the set
$Y$ there is an open hole. Countable union of nowhere dense sets is called the
set of first category. The second usual way is to define thin sets as the sets
of first category. A residual set (a “thick" set) is the complement of a set
of the first category. If a set is not meagre it is said to be of the second
category. The Baire classification is non-trivial in the so-called Baire
spaces, where every intersection of a countable collection of open dense sets
is also dense. Complete metric spaces and, in particular, Banach spaces are
Baire spaces. Therefore, for Banach spaces of functions, the common definition
of negligible set is “set of first Baire category”. Such famous results as
transversality theorem in differential topology Golub1974  or Pugh closing
lemma Pugh1967  and Kupka-Smale theorem Palis1982  in differential dynamics.

Despite of these great successes, it is also widely recognized that the Bair
category approach to generic properties requires at least great care. Here are
some examples of correct but useless statements about “generic" properties of
function: Almost every continuous function is not differentiable; Almost every
$C^{1}$-function is not convex. Their meaning for applications is most
probably this: the genericity used above for continuous functions or for
$C^{1}$-function is irrelevant to the subject.

Contradictions between the measure-based and category-based definitions of
negligible sets are well-known even in dimension one: even the real line $R$
can be divided into two sets, one of which has zero measure, the other is of
first category Oxtoby2013 . Genericity in the sense of measure and genericity
in the sense of category differ significantly in the applications where both
concepts can be used.

The conflict between the two main views on genericity and negligibility
stimulated efforts to invent new and stronger approaches. The formal
requirements to new definitions are:

  * •

A union of countable family of thin sets should be thin.

  * •

Any subset of a thin set should be thin.

  * •

The whole space is not thin.

Of course, if we take care not to throw the baby out with the bath water then
in $\mathbb{R}^{n}$, where both classical definition are applicable, we expect
that thin sets should be of first category and have zero measure, both. It was
not clear a priori whether such a theory is possible with proof nontrivial and
important generic properties. It turned out that it is possible. To
substantiate the effectiveness of evolutionary optimization, a theory of
completely negligible sets in Banach spaces was developed. Gorban1984 ;
Gorban2007 .

Let $Q$ be a real Banach space. Consider compact subsets in $Q$ parametrized
by points of a compact space $K$. It can be presented as a Banach space
$C(K,Q)$ of continuous maps $K\rightarrow Q$ in the maximum norm.

A set $Y\subset Q$ is completely thin, if for any compact space $K$ the set of
continuous maps $\Psi:K\rightarrow Q$ with non-empty intersection $\Psi(K)\cap
Y\neq\emptyset$ is set of first Bair category.

The union of a finite or countable number of completely thin sets is
completely thin. Any subset of a completely thin point is completely thin,
while the whole $Q$ is not. A set $Y$ in the Banach space $Q$ is completely
thin, if for any compact set $K$ in $Q$ and arbitrary positive $\varepsilon>0$
there exists a vector $q\in Q$, such that $\|q\|<\varepsilon$ and $K+q$ does
not intersect $Y$: $(K+q)\cap Y=\emptyset$. All compact sets in infinite-
dimensional Banach spaces and closed linear subspaces with infinite
codimension are completely thin.

Only empty set is completely thin in a finite-dimensional space
$\mathbb{R}^{n}$.

Examples below demonstrate that almost all continuous functions have very
natural properties: the set of zeros is nowhere dense, and the (global)
maximizer is unique. Below the wording “almost always" means: the set of
exclusions is completely thin.

[Gorban1984 ; Gorban2007 ]Let $X$ have no isolated points. Then

  * •

Almost always a function $f\in C(X)$ has nowhere dense set of zeros $\\{x\in
X\,|\,f(x)=0\\}$ (the set of exclusions is completely thin in $C(X)$).

  * •

Almost always a function $f\in C(X)$ has only one point of global maximum.

The following proposition is a tool for proof that some typical properties of
functions hold almost always for all functions from a generic compact set.

[Gorban1984 ; Gorban2007 ]If a set $Y$ in the Banach space $Q$ is completely
thin, then for any compact metric space $K$ the set of continuous maps
$\Psi:K\rightarrow Q$ with non-empty intersection $\Psi(K)\cap Y\neq\emptyset$
is completely thin in the Banach space $C(K,Q)$.

[Gorban1984 ; Gorban2007 ] Let $X$ have no isolated points. Then for any
compact space $K$ and almost every continuous map $\Psi:K\rightarrow C(X)$ all
functions $f\in\Psi(K)$ have nowhere dense sets of zeros (the set of
exclusions is completely thin in $C(K,C(X))$).

In other words, in almost every compact family of continuous functions all the
functions have nowhere dense sets of zeros.

Qualitatively, the concept of a completely thin set was introduced as a tool
for identifying typical properties of infinite-dimensional objects, the
violation of which is unlikely (‘improbable’) in any reasonable sense.

##  Appendix B Flowchart of multi-corrector operation

In Sec. 4, we introduced multi-corrector of AI systems. The basic scheme of
this device is presented in Fig. 6. It includes several elementary correctors
(see Fig. 2) and a dispatcher. A cluster of errors is owned by each elementary
corrector. An elementary corrector evaluates the risk of errors from its own
cluster for an arbitrary operation situation and takes the decision: to
correct or not to correct the legacy AI decision for this situation. For any
situation, the dispatcher selects the most appropriate elementary corrector to
make a decision about correction. To find a suitable corrector, it uses a
cluster error model. When new errors are found, the cluster model changes.
More detailed presentation of multi-corrector operation is given by the
following flowcharts.

Figure 10:  Notations used in the flowcharts. All flowcharts use a unified set
of blocks: blocks in the form of parallelograms display data, rectangular
blocks display procedures, and blocks in the form of rhombuses display the
branching points of processes (algorithms) or decision points. The arrows
reflect the transfer of data and control.

Flowcharts and blocks are numbered. The flowchart number is mentioned at the
top of the drawings. If a block is present in different flowcharts, then it
carries the number assigned to it in the top-level flowchart. The relations
between different flowcharts are presented in Fig. 11.

Figure 11:  The tree of flowcharts: 10 – Operation of the modified AI system
(Fig. 12); 12 – Operation of the legacy AI system (Fig. 13), 14 – Operation of
the correction system (Fig. 14); 143 – Single corrector operation (Fig. 15);
141 – The work of the dispatcher (Fig. 16); 1415 – Online modification of the
cluster model (Fig. 17);  Figure 12:  Operation of the modified AI system
(10). Input signals (11) are fed to the input of the AI system (12), which at
the output gives out the complete vector of the signal (13) that can be used
for correction. The complete signal vector (13) is fed to the input of the
correction system (14). The correction system (14) calculates the correction
of the output signals (15).  Figure 13:  Operation of the legacy AI system
(12). Input signals (11) are fed to the input of the AI system. The AI system
generates vectors of internal signals (123) and output signals (122). Input
signals (11), internal signals (123), and output signals (122) form the
complete signal vector (13).  Figure 14:  Operation of the correction system
(14). The complete vector of signals (13) is fed to the dispatcher input
(141). The dispatcher (141) selects from the correctors the one that most
closely matches the situation (142). The selected corrector (142) and the
complete signal vector (13) are used to correct the signals (13). The computed
corrected outputs (15) are returned.  Figure 15:  Single corrector operation
(143). The complete vector of signals (13) is used to decide whether a
correction is needed (1431). If it is necessary, then correction (1435) is
performed, and the resulting vector of output signals (15) is sent to the
output. If there is no need for correction, then the vector of output signals
is extracted (1434) from the complete vector of signals (13), and the
resulting vector of output signals (15) is transmitted to the output.  Figure
16: The work of the dispatcher (141). If the error flag (1412) is detected
(1414), then the current cluster model (1411) and the complete signal vector
(13) are used to modify the cluster model (1415) online. The modified cluster
model becomes the current one (1411). If the error flag (1412) is not detected
(1414), then the current cluster model (1411) is selected (1413) for use
(1411). Based on the cluster model (1411) and the complete signal vector (13),
the most suitable cluster (1416) is selected. Then, the corrector (142)
corresponding to this cluster is selected (1417).  Figure 17:  Online
modification of the cluster model (1415). Based on the current cluster model
(1411) and the complete signal vector (13), the most suitable cluster (14151)
is selected (1416). Online modification of the rule for determining the most
suitable cluster (14152) is performed. After setting up the new cluster model
(1411), the classifier for this corrector to make a decision about the need
for correction is explicitly modified (14153). The modified corrector (142)
together with the new cluster model (1411) forms an updated version of the
correction system (14).

References

## References

  * (1) Vinyals, O.; Blundell, C.; Lillicrap, T.; Kavukcuoglu, K.; Wierstra, D. Matching networks for one shot learning. NIPS’16: Proceedings of the 30th International Conference on Neural Information Processing Systems, December 2016, Pages 3637–3645. https://dl.acm.org/doi/abs/10.5555/3157382.3157504
  * (2) Ravi, S.; Larochelle, H. Optimization as a model for few-shot learning. Proceedings of the International Conference on Learning Representations, 2017. https://openreview.net/pdf?id=rJY0-Kcll
  * (3) Wang, Y.; Yao, Q.; Kwok, J. T.; Ni, L. M. Generalizing from a few examples: A survey on few-shot learning. ACM Computing Surveys (CSUR) 2020, 53(3), 1–34. https://doi.org/10.1145/3386252
  * (4) Snell, J.; Swersky, K.; Zemel, R. Prototypical networks for few-shot learning. In Proceedings of the 31st International Conference on Neural Information Processing Systems, 2017, 4080–4090. https://proceedings.neurips.cc/paper/2017/file/cb8da6767461f2812ae4290eac7cbc42-Paper.pdf
  * (5) Sung, F.; Yang, Y.; Zhang, L.; Xiang, T.; Torr, P.H.; Hospedales, T. M. Learning to compare: Relation network for few-shot learning. In Proceedings of the IEEE conference on computer vision and pattern recognition, 2018, 1199–1208). https://openaccess.thecvf.com/content_cvpr_2018/html/Sung_Learning_to_Compare_CVPR_2018_paper.html
  * Gorban and Tyukin (2018) Gorban, A.N.; Tyukin, I.Y.  Blessing of dimensionality: mathematical foundations of the statistical physics of data.  Phil. Trans. R. Soc. A 2018, 376, 20170237, https://doi.org/10.1098/rsta.2017.0237.
  * (7) Tyukin, I.Y.; Gorban, A.N.; Alkhudaydi, M.H.; Zhou, Q. Demystification of few-shot and one-shot learning. arXiv preprint arXiv:2104.12174 2021. https://arxiv.org/abs/2104.12174
  * (8) Kainen, P.C.  Utilizing geometric anomalies of high dimension: when complexity makes computation easier. In Computer-Intensive Methods in Control and Signal Processing: The Curse of Dimensionality; Warwick, K., Kárný, M., Eds.; Springer: New York, NY, USA, 1997; 283–294, https://doi.org/10.1007/978-1-4612-1996-5_18.
  * Donoho (2000) Donoho, D.L.  _High-Dimensional Data Analysis: The Curses and Blessings of Dimensionality_ ;  Invited lecture at Mathematical Challenges of the 21st Century, AMS National Meeting, Los Angeles, CA, USA, August 6-12, 2000; CiteSeerX http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.329.3392.
  * Anderson et al. (2014) Anderson, J.; Belkin, M.; Goyal, N.; Rademacher, L.; Voss, J.  The More, the Merrier: the Blessing of Dimensionality for Learning Large Gaussian Mixtures.  In Proceedings of The 27th Conference on Learning Theory, Barcelona, Spain, 13–15 June 2014; Balcan, M.F.; Feldman, V.; Szepesvári, C., Eds.; PMLR: Barcelona, Spain, 2014; Volume 35, 1135–1164. http://proceedings.mlr.press/v35/anderson14.pdf
  * Gorban et al. (2016) Gorban, A.N.; Tyukin, I.Y.; Romanenko, I.  The blessing of dimensionality: Separation theorems in the thermodynamic limit.  IFAC-PapersOnLine 2016, 49, 64–69, https://doi.org/10.1016/j.ifacol.2016.10.755.
  * (12) Gorban, A.N.; Tyukin, I.Y. Stochastic separation theorems. Neural Networks 2017, 94, 255–259. https://doi.org/10.1016/j.neunet.2017.07.014
  * Gorban et al. (2018) Gorban, A.N.; Golubkov, A.; Grechuk, B.; Mirkes, E.M.; Tyukin I.Y. Correction of AI systems by linear discriminants: Probabilistic foundations. Information Sciences 2018, 466, 303–322. https://doi.org/10.1016/j.ins.2018.07.040
  * Gorban et al. (2019) Gorban, A.N.; Makarov, V.A.; Tyukin, I.Y.  The unreasonable effectiveness of small neural ensembles in high-dimensional brain.  Phys. Life Rev. 2019, 29, 55–88, https://doi.org/10.1016/j.plrev.2018.09.005.
  * (15) Grechuk, B.; Gorban, A.N.; Tyukin, I.Y. General stochastic separation theorems with optimal bounds. Neural Networks 2021, 138, 33–56. https://doi.org/10.1016/j.neunet.2021.01.034
  * (16) Flury B. Principal points. Biometrika 1990, 77, 33–41. https://doi.org/10.1093/biomet/77.1.33
  * (17) Gorban, A.N.; Zinovyev, A.Y. Principal graphs and manifolds. In Handbook of research on machine learning applications and trends: algorithms, methods, and techniques, Olivas, E.S.; Guererro, J.D.M.; Sober, M.M.; Benedito, J.R.M.; Lopes, A.J.S., Eds. 2010, 28–59. IGI Global. https://doi.org/10.4018/978-1-60566-766-9.ch002
  * (18) Tyukin, I.Y.; Gorban, A.N.; Grechuk, B.; Green, S. Kernel Stochastic Separation Theorems and Separability Characterizations of Kernel Classifiers. In 2019 International Joint Conference on Neural Networks (IJCNN) 2019, 1–6, IEEE. https://doi.org/10.1109/IJCNN.2019.8852278
  * (19) Kreinovich, V.; Kosheleva, O. Limit Theorems as Blessing of Dimensionality: Neural-Oriented Overview. Entropy 2021, 23(5), 501. https://doi.org/10.3390/e23050501
  * Gorban et al. (2008) Gorban, A.N.; Kégl, B.; Wunsch, D.; Zinovyev, A. (Eds.)  Principal Manifolds for Data Visualisation and Dimension Reduction; Springer: Berlin/Heidelberg, Germany, 2008. https://doi.org/10.1007/978-3-540-73750-6.
  * (21) Schölkopf, B. Nonlinear Component Analysis as a Kernel Eigenvalue Problem. Neural Computation 1998, 10(5), 1299–1319. https://doi.org/10.1162/089976698300017467.
  * (22) Gorban, A.N.; Zinovyev, A. Principal manifolds and graphs in practice: from molecular biology to dynamical systems. International Journal of Neural Systems 2010, 20(03), 219–232. https://doi.org/10.1142/S0129065710002383
  * (23) Kramer, M.A. Nonlinear principal component analysis using autoassociative neural networks. AIChE Journal 199137(2), 233–243. https://doi.org/10.1002/aic.690370209
  * (24) Hinton G.E.; Salakhutdinov, R.R. Reducing the dimensionality of data with neural networks. Science 2006, 28313(5786), 504–507. https://doi.org/10.1126/science.1127647
  * (25) Gorban, A.N.; Makarov, V.A.; Tyukin, I.Y. High-Dimensional Brain in a High-Dimensional World: Blessing of Dimensionality. Entropy 2020, 22, 82, https://doi.org/10.3390/e22010082.
  * Giannopoulos and Milman (2000) Giannopoulos, A.A.; Milman, V.D.  Concentration property on probability spaces.  Adv. Math. 2000, 156, 77–106, https://doi.org/10.1006/aima.2000.1949.
  * Gromov (2003) Gromov, M.  Isoperimetry of waists and concentration of maps.  Geom. Funct. Anal. 2003, 13, 178–215, https://doi.org/10.1007/s00039-009-0703-1.
  * Ledoux (2005) Ledoux, M.  The Concentration of Measure Phenomenon; Number 89 in Mathematical Surveys & Monographs; AMS: Providence, RI, USA, 2005.
  * Vershynin (2018) Vershynin, R.  High-Dimensional Probability: An Introduction with Applications in Data Science; Cambridge Series in Statistical and Probabilistic Mathematics; Cambridge University Press: Cambridge, UK, 2018.
  * (30) Gartner Hype Cycle for Artificial Intelligence, 2019. https://www.gartner.com/smarterwithgartner/top-trends-on-the-gartner-hype-cycle-for-artificial-intelligence-2019/.
  * (31) Gartner Hype Cycle for Emerging Technologies, 2020. https://www.gartner.com/en/newsroom/press-releases/2020-08-18-gartner-identifies-five-emerging-trends-that-will-drive-technology-innovation-for-the-next-decade.
  * (32) Gorban, A.N.; Grechuk, B.; Tyukin, I.Y. Augmented Artificial Intelligence: a Conceptual Framework, arXiv preprint 2018. https://arxiv.org/abs/1802.02172.
  * Kainen & Kůrková (1993) Kainen, P.; Kůrková, V. Quasiorthogonal dimension of Euclidian spaces. Appl. Math. Lett. 1993, 6, 7–10. https://doi.org/10.1016/0893-9659(93)90023-G.
  * Kainen & Kůrková (2020) Kainen, P.; Kůrková, V. Quasiorthogonal dimension. In Kosheleva, O., Shary, S.P., Xiang, G., Zapatrin, R. (Eds.). Beyond Traditional Probabilistic Data Processing Techniques: Interval, Fuzzy etc. Methods and Their Applications. Springer, Cham, 2020, 615–629. https://doi.org/10.1007/978-3-030-31041-7_35.
  * Gorban et al. (2016a) Gorban, A.N.; Tyukin, I.; Prokhorov, D.; Sofeikov, K. Approximation with random bases: Pro et contra. Information Sciences 2016 364–365, 129–145. https://doi.org/10.1016/j.ins.2015.09.021
  * (36) Camastra F. Data dimensionality estimation methods: a survey. Pattern recognition 2003, 36(12), 2945–2954. https://doi.org/10.1016/S0031-3203(03)00176-6
  * (37) Bac J.; Zinovyev A. Lizard brain: Tackling locally low-dimensional yet globally complex organization of multi-dimensional datasets. Frontiers in neurorobotics 2020 13, 110. https://doi.org/10.3389/fnbot.2019.00110
  * (38) Moczko, E.; Mirkes, E.M.; Cáceres, C.; Gorban, A.N.; Piletsky, S. Fluorescence-based assay as a new screening tool for toxic chemicals. Scientific reports 2016, 6(1), 33922. https://doi.org/10.1038/srep33922
  * (39) Chen, S. H.; Pollino, C. A. Good practice in Bayesian network modelling. Environmental Modelling & Software 2012, 37, 134–145. https://doi.org/10.1016/j.envsoft.2012.03.012
  * (40) Cobb, B.R.; Rumí, R.; Salmerón, A. Bayesian Network Models with Discrete and Continuous Variables. In: Lucas P., Gámez J.A., Salmerón A. (eds) Advances in Probabilistic Graphical Models. Studies in Fuzziness and Soft Computing, vol 213. Springer, Berlin, Heidelberg, 2007. https://doi.org/10.1007/978-3-540-68996-6_4
  * (41) Rezende, D.J.; Mohamed, S.; Wierstra, D. Stochastic backpropagation and approximate inference in deep generative models. In Proceedings of the 31st International Conference on Machine Learning 2014 Jun 22-24 Eds. E.P. Xing and T. Jebara (pp. 1278–1286). PMLR. http://proceedings.mlr.press/v32/rezende14.html
  * (42) Hinton G.E. A Practical Guide to Training Restricted Boltzmann Machines. In: Montavon G., Orr G.B., Müller KR. (eds) Neural Networks: Tricks of the Trade. Lecture Notes in Computer Science, vol 7700. Springer, Berlin, Heidelberg, 2012, 599–619. https://doi.org/10.1007/978-3-642-35289-8_32
  * (43) Streiner, D.L.; Norman, G.R. Correction for multiple testing: is there a resolution?. Chest 2011140(1), 16–18. https://doi.org/10.1378/chest.11-0523
  * (44) Noble, W.S. How does multiple testing correction work?. Nat Biotechnol 2009, 27, 1135–1137. https://doi.org/10.1038/nbt1209-1135
  * Jolliffe (1993) Jolliffe, I. Principal Component Analysis; Springer: Berlin/Heidelberg, 1993.
  * (46) Sompairac, N.; Nazarov, P.V.; Czerwinska, U.; Cantini, L.; Biton, A.; Molkenov, A.; Zhumadilov, Z.; Barillot, E.; Radvanyi, F.; Gorban, A.; Kairov, U. Independent component analysis for unraveling the complexity of cancer omics datasets. International Journal of molecular sciences 2019, 20(18), 4414. http://dx.doi.org/10.3390/ijms20184414
  * (47) Hicks, S.C.; Townes, F.W.; Teng, M.; Irizarry, R.A.. Missing data and technical variability in single-cell RNA-sequencing experiments. Biostatistics 2018, 19(4), 562–578. https://doi.org/10.1093/biostatistics/kxx053
  * (48) Krumm, N.; Sudmant, P.H.; Ko, A.; O’Roak, B.J.; Malig, M.; Coe, B.P.; Quinlan, A.R.; Nickerson, D.A.; Eichler, E.E., Copy number variation detection and genotyping from exome sequence data. Genome research 22(8), 1525–1532. https://doi.org/10.1101/gr.138115.112
  * (49) Koren, Y.; Carmel, L. Robust linear dimensionality reduction. IEEE Trans Visual Comput Graph. 2004 10(4), 459–470. https://doi.org/10.1109/TVCG.2004.17
  * (50) Mirkes, E.M., Gorban, A.N., Zinoviev, A. Supervised PCA. https://github.com/Mirkes/SupervisedPCA. 2016.
  * (51) Gorban, A.N.; Mirkes, E.M.; Tukin I.Y. How deep should be the depth of convolutional neural networks: a backyard dog case study. Cognitive Computation 2020, 12, 388–397. https://doi.org/10.1007/s12559-019-09667-7
  * (52) Song, Y.; Nie, F.; Zhang, C.; Xiang, S. A unified framework for semi-supervised dimensionality reduction. Pattern Recognition 2008, 41(9), 2789–2799. https://doi.org/10.1016/j.patcog.2008.01.001.
  * (53) Cangelosi, R.; Goriely, A. Component retention in principal component analysis with application to cDNA microarray data. Biol. Direct 2007, 2(1), 2. doi:10.1186/1745-6150-2-2
  * (54) Ben-David, S.; Blitzer, J.; Crammer, K.; Kulesza, A.; Pereira, F.; Vaughan, J.W. A theory of learning from different domains. Mach. Learn. 2010, 79, 151–175. https://doi.org/10.1007/s10994-009-5152-4
  * (55) Sun, S.: Shi, H.; Wu, Y. A survey of multi-source domain adaptation. Information Fusion 2015, 24, 84–92. https://doi.org/10.1016/j.inffus.2014.12.003
  * (56) Saito, K.; Watanabe, K., Ushiku, Y.; Harada, T. Maximum classifier discrepancy for unsupervised domain adaptation. In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition 2018, 3723–3732. https://doi.org/10.1109/CVPR.2018.00392
  * (57) Ganin, Y.; Ustinova, E.; Ajakan, H.; Germain, P.; Larochelle, H.; Laviolette, F.; Marchand, M.; Lempitsky, V. Domain-adversarial training of neural networks. The Journal of Machine Learning Research 2016, 17(59), 2096–2030. http://jmlr.org/papers/v17/15-239.html
  * (58) Matasci, G.; Volpi, M.; Tuia, D.; Kanevski, M. Transfer component analysis for domain adaptation in image classification. In Image and Signal Processing for Remote Sensing XVII (Vol. 8180, p. 81800F). International Society for Optics and Photonics, 2011. https://doi.org/10.1117/12.898229
  * (59) Pestov, V. Is the k-NN classifier in high dimensions affected by the curse of dimensionality? Comput. Math. Appl. 2013, 65, 1427–1437. https://doi.org/10.1016/j.camwa.2012.09.011
  * (60) Mirkes, E.M.; Allohibi, J.; Gorban, A.N. Fractional Norms and Quasinorms Do Not Help to Overcome the Curse of Dimensionality. Entropy 2020, 22, 1105. https://doi.org/10.3390/e22101105.
  * (61) Zadeh, L.A. Towards a theory of fuzzy information granulation and its centrality in human reasoning and fuzzy logic. Fuzzy Sets Syst. 1997, 19, 111–127. https://doi.org/10.1016/S0165-0114(97)00077-8
  * (62) Pedrycz, W.; Skowron, A.; Kreinovich, V. (Eds) Handbook of granular computing. John Wiley & Sons; 2008.
  * (63) Kainen, P.C. Replacing points by compacta in neural network approximation. Journal of the Franklin Institute 2004, 341(4), 391–399. https://doi.org/10.1016/j.jfranklin.2004.03.001
  * Guédon and Milman (2011) Guédon, O.; Milman, E. Interpolating thin-shell and sharp large-deviation estimates for lsotropic log-concave measures. Geometric and Functional Analysis 2011, 21(5), 1043–1068. https://doi.org/10.1007/s00039-011-0136-5
  * Lévy (1951) Lévy, P.  Problèmes Concrets D’analyse Fonctionnelle. Gauthier-Villars: Paris, France, 1951.
  * (66) Khinchin, A.Y. Mathematical Foundations of Statistical Mechanics. New York: Courier Corporation, 1949. (English translation from the Russian edition, Moscow – Leningrad, 1943.)
  * (67) Thompson, C.J. Mathematical statistical mechanics. Princeton University Press; 2015.
  * (68) Kolmogorov, A.N., Foundations of the theory of probability. Second English Edition. Courier Dover Publications; 2018.
  * (69) Liu, L.; Shao, L.; Li, X. Evolutionary compact embedding for large-scale image classification. Information Sciences 2015, 316, 567–581. https://doi.org/10.1016/j.ins.2014.06.030
  * (70) Vemulapalli, R.; Agarwala, A. A compact embedding for facial expression similarity. In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition 2019, 5683–5692. https://openaccess.thecvf.com/content_CVPR_2019/html/Vemulapalli_A_Compact_Embedding_for_Facial_Expression_Similarity_CVPR_2019_paper.html
  * (71) Bhattarai, B.; Liu, H.; Huang, H.H. Ceci: Compact embedding cluster index for scalable subgraph matching. In Proceedings of the 2019 International Conference on Management of Data 2019 Jun 25 (pp. 1447-1462). https://doi.org/10.1145/3299869.3300086
  * Tyukin et al. (2020) Tyukin, I.Y.; Higham, D.J.; Gorban, A. N. On adversarial examples and stealth attacks in artificial intelligence systems. In Proc. 2020 International Joint Conference on Neural Networks (IJCNN), Glasgow, United Kingdom, 2020 (pp. 1–6), IEEE, 2020. https://doi.org/10.1109/IJCNN48605.2020.9207472.
  * (73) Colbrook, M.J.; Antun. V.; Hansen, A.C. Can stable and accurate neural networks be computed?–On the barriers of deep learning and Smale’s 18th problem. arXiv preprint arXiv:2101.08286 2021.
  * (74) Rudin, W. Functional Analysis. New York, NY: McGraw-Hill Science/Engineering/Math, 1991.
  * (75) Xu, R; Wunsch, D. Clustering. John Wiley & Sons; 2008.
  * (76) Tyukin, I.Y., Gorban, A.N., McEwan, A.A., Meshkinfamfard, S.; Tang, L. Blessing of dimensionality at the edge and geometry of few-shot learning. Information Sciences 2021, 564, 124–143. https://doi.org/10.1016/j.ins.2021.01.022
  * (77) Tao, C.W. Unsupervised fuzzy clustering with multi-center clusters. Fuzzy Sets and Systems 2002, 128(3), 305–322. https://doi.org/10.1016/S0165-0114(01)00191-9
  * (78) Krizhevsky, A. Learning multiple layers of features from tiny images. Technical Report, University of Toronto, 2009; CiteSeerX https://citeseerx.ist.psu.edu/viewdoc/versions?doi=10.1.1.222.9220
  * (79) Krizhevsky, A. CIFAR 10 Dataset, University of Toronto, 2009. https://www.cs.toronto.edu/~kriz/cifar.html
  * (80) Golubitsky, M., Guillemin, V. Stable Mappings and Their Singularities. Springer-Verlag. 1974.
  * (81) Pugh, C. The closing lemma. Amer. J. Math. 1967, 89, 956–1009. https://doi.org/10.2307/2373413
  * (82) Palis, J.; de Melo W. The Kupka-Smale Theorem. In: Geometric Theory of Dynamical Systems. Springer, New York, NY, 1982. https://doi.org/10.1007/978-1-4612-5703-5_3
  * (83) Oxtoby, J.C. Measure and category: A survey of the analogies between topological and measure spaces. Springer; NY, 2013.
  * (84) Gorban, A.N. Equilibrium encircling. Equations of chemical kinetics and their thermodynamic analysis. Nauka, Novosibirsk, 1984.
  * (85) Gorban, A.N. Selection Theorem for Systems with Inheritance, Math. Model. Nat. Phenom. 2007, 2(4), 1–45. https://doi.org/10.1051/mmnp:2008024
  * (86) Feynman, R.P.; Leighton, R.B.; Sands, M. The Feynman Lectures on Physics, Addison Wesley, 2005, Volume I, Chapter 1, Introduction.


Discrete Bessel and Mathieu functions

Kenan Uriostegui111Posgrado en Ciencias Físicas, Universidad Nacional Autónoma
de México and Kurt Bernardo Wolf

Instituto de Ciencias Físicas

Universidad Nacional Autónoma de México

Av. Universidad s/n, Cuernavaca, Morelos 62251, México

###### Abstract

The two-dimensional Helmholtz equation separates in elliptic coordinates based
on two distinct foci, a limit case of which includes polar coordinate systems
when the two foci coalesce. This equation is invariant under the Euclidean
group of translations and orthogonal transformations; we replace the latter by
the discrete dihedral group of $N$ discrete rotations and reflections. The
separation of variables in polar and elliptic coordinates is then used to
define discrete Bessel and Mathieu functions, as approximants to the well-
known continuous Bessel and Mathieu functions, as $N$-point Fourier transforms
approximate the Fourier transform over the circle, with integrals replaced by
finite sums. We find that these ‘discrete’ functions approximate the numerical
values of their continuous counterparts very closely and preserve some key
special function relations.

##  1 Introduction

The role of the Euclidean group of translations, reflections and rotations in
the determination of the coordinate systems that separate the solutions of the
two-dimensional Helmholtz equation is well known from the work by Willard
Miller Jr. [1, Ch. 1]. This symmetry accounts for their separability in four
coordinate systems: Cartesian, polar, parabolic and elliptic. Only the
elliptic system is generic; when the two foci coalesce, this system becomes
the polar one with angular and radial coordinates; when one focus departs to
infinity the system becomes parabolic; and when both foci do, it becomes
Cartesian.

The polar decomposition was used by Biagetti et al. [2] to first introduce a
discrete version of Bessel functions based on an expansion of plane waves into
a finite number of polar components —that was not quite complete. This was
properly completed in Ref. [3], defining discrete Bessel functions
$B_{n}^{{\scriptscriptstyle N}}(\rho)$, which approximate the usual continuous
Bessel functions $J_{n}(\rho)$ by replacing Fourier series over a circle
${\cal S}^{1}$ by the finite Fourier transform on $N$ equidistant points on
that circle,

$\theta_{m}=2\pi m/N,\qquad m\in\\{0,1,\ldots,N{-}1\\}=:{\cal
S}^{1}_{{\scriptscriptstyle(N)}},$ (1)

where $m$ is counted modulo $N$. It was found that these discrete functions
approximated very closely (of the order $10^{-16}$) the corresponding
continuous ones over a region, roughly $0\leq n+\rho<N$.

Several authors have introduced functions that approximate the well-known
continuous Bessel functions $J_{n}(\rho)$, for the purpose of reducing
computation time, or to provide new classes of solutions to difference
equations that will share some of their salient properties [4, 5, 6]. Our
approach follows the well known approximation afforded by the $N$-point finite
Fourier transform to the integral Fourier transform over the circle. This is
done for polar and elliptic coordinates, and introduces both ‘discrete’ Bessel
and Mathieu functions. These functions, we should emphasize, differ from those
proposed in the works cited above, which are also distinct in definition and
purpose among themselves. By construction it will follow that under
$N\to\infty$, these discrete functions become the continuous ones, although
this limit requires further mathematical precision, as it may involve Gibbs-
type oscillation phenomena that we cannot address here.

In Sect. 2 we present this discretization method and a resumé of the results
in Ref. [3] for Bessel functions, to note that the discrete functions thus
defined approximate the continuous ones remarkably well. In Sect. 3 of the
present paper we apply the strategy of replacing harmonic analysis on ${\cal
S}^{1}$ by ${\cal S}^{1}_{{\scriptscriptstyle(N)}}$ to define discrete
approximants to the Mathieu functions of first and second kind in the elliptic
coordinate system. All relations are backed by numerical verification. In the
concluding Sect. 4 we provide some further connections and preliminary
conclusions.

##  2 Continuous and discrete Bessel functions

The Helmholtz equation for wavefields $f(x,y)$ of (fixed) real wavenumber
$\kappa\in{\cal R}$, is

$(\partial_{x}^{2}+\partial_{y}^{2}+\kappa^{2})f(x,y)=0,$ (2)

with $\partial_{z}\equiv\partial/\partial_{z}$ and $(x,y)\in{\cal R}^{2}$. In
this section we follow the well known case of polar coordinates,

$x=r\cos\theta,\quad y=r\sin\theta,\qquad
r\in[0,\infty),\quad\theta\in(-\pi,\pi]={\cal S}^{1}.$ (3)

A key assumption is a Hilbert space structure for the solutions $f(x,y)$ by
which one can write them as the two-dimensional Fourier transform,

$f(x,y)=\frac{1}{2\pi}\int\\!\\!\\!\int_{{\cal R}^{2}}{\rm d}\kappa_{x}\,{\rm
d}\kappa_{y}\,\exp{\rm
i}(x\kappa_{x}+y\kappa_{y})\,\widetilde{f}(\kappa_{x},\kappa_{y}).$ (4)

The Helmholtz equation (2) is then correspondingly transformed to a conjugate
space $(\kappa_{x},\kappa_{y})\in{\cal R}^{2}$ where it reads
$(\kappa^{2}-\kappa_{x}^{2}-\kappa_{y}^{2})\widetilde{f}(\kappa_{x},\kappa_{y})=0$,
which we can also refer to polar coordinates $\kappa_{x}=\kappa\cos\phi$,
$\kappa_{y}=\kappa\sin\phi$, with the surface element ${\rm d}\kappa_{x}\,{\rm
d}\kappa_{y}=\kappa\,{\rm d}\kappa\,{\rm d}\phi$. The solutions to the
Fourier-transformed Helmholtz equation are thus reduced by a Dirac
$\delta$-distributions in the radius [1, Ch. 1], as
$\widetilde{f}(\kappa_{x},\kappa_{y})=\sqrt{2\pi}\kappa^{-1}\delta(\kappa-\widetilde{\kappa})\,{f}_{\\!\circ}(\phi)$,
with ${f}_{\\!\circ}(\phi)$ a function on the $\phi$-circle ${\cal S}^{1}$ of
radius $\widetilde{\kappa}$, that we write again $\kappa$, understanding that
it is the fixed wavenumber. The Helmholtz solutions (4) thus acquire the
single-integral form

$f(x,y)=\frac{1}{\surd 2\pi}\int_{{\cal S}^{1}}{\rm d}\phi\,\exp{\rm
i}\kappa(x\cos\phi+y\sin\phi)\,{f}_{\\!\circ}(\phi),$ (5)

with the Hilbert space structure based on the inner product of functions
$f^{{\scriptscriptstyle(1)}}_{\\!\circ}(\phi)$ and
$f^{{\scriptscriptstyle(2)}}_{\\!\circ}(\phi)$ on the circle,

$(f^{{\scriptscriptstyle(1)}}_{\\!\circ},f^{{\scriptscriptstyle(2)}}_{\\!\circ})_{\circ}:=\int_{{\cal
S}^{1}}{\rm
d}\phi\,f^{{\scriptscriptstyle(1)}}_{\\!\circ}(\phi)^{*}f^{{\scriptscriptstyle(2)}}_{\\!\circ}(\phi).$
(6)

It is here that we reduce the continuous circle Fourier transform to the
$N$-point discrete Fourier transform, from ${{\cal S}^{1}}$ to ${\cal
S}^{1}_{{\scriptscriptstyle(N)}}$, replacing integrals by summations and the
continuous variable $\phi\in{\cal S}^{1}$ with $\phi_{m}\in{\cal
S}^{1}_{{\scriptscriptstyle(N)}}$, as

$\int_{{\cal S}^{1}}{\rm
d}\phi\,{F}_{\\!\circ}(\phi)\leftrightarrow\sum_{m\in{\cal
S}^{1}_{{\scriptscriptstyle(N)}}}{F}(\phi_{m}),\qquad{2\pi\leftrightarrow
N,\atop\phi_{m}=2\pi m/N,}$ (7)

for $m\in\\{0,1,\ldots,N{-}1\\}$ counted modulo $N$; the set of $N$ discrete
angles $\phi_{m}$ are thus equidistant by $2\pi/N$. The functions
$f(\phi_{m})\equiv f_{m}$ can be interpreted as sample points of a continuous
function, or as the index for the list of components of an $N$-cyclic vector.
In either case, the inner product of two discrete functions
$f^{{\scriptscriptstyle(1)}}_{n}$ and $f^{{\scriptscriptstyle(2)}}_{n}$ is
naturally

$(f^{{\scriptscriptstyle(1)}},f^{{\scriptscriptstyle(2)}})_{{\scriptscriptstyle(N)}}:=\sum_{n=0}^{N-1}f_{n}^{{{\scriptscriptstyle(1)}}*}\,f^{{\scriptscriptstyle(2)}}_{n},$
(8)

and it is clear that the $N\to\infty$ limit will lead back from the discrete
to the continuum, with the approximations and limits familiar from Fourier
theory.

The Helmholtz equation (2) in polar coordinates, multiplied by $r^{2}$,

$(r^{2}\partial^{2}_{r}+r\partial_{r}+\partial^{2}_{\phi}+\kappa^{2})f(r,\phi)=0,$
(9)

shows that solutions can be factored into a function of the radius times a
function of the angle as $f(r,\phi)=R(r)\,\Phi(\phi)$, while (5) implies that
solutions $\Phi(\phi)$ for the angular factor will determine a corresponding
radial factor $R(r)$. An orthonormal and complete set of eigenfunctions of
$\partial^{2}_{\phi}$ over the circle $\phi\in{\cal S}^{1}$ is the set of
phases $\Phi_{n}(\phi):=(2\pi)^{-1/2}\exp({\rm i}n\phi)$, with integer
$n\in\\{0,\pm 1,\ldots\\}$, and inner products
$(\Phi_{n},\Phi_{n^{\prime}})_{\circ}=\delta_{n,n^{\prime}}$. When the domain
of these functions is restricted from $\phi\in{\cal S}^{1}$ to
$\phi_{m}\in{\cal S}^{1}_{{\scriptscriptstyle(N)}}$ as in (1), we retain the
subset of $N$ functions on the $N$ points in ${\cal
S}^{1}_{{\scriptscriptstyle(N)}}$, given by

$\Phi^{{\scriptscriptstyle(N)}}_{n}(\phi_{m}):=\frac{1}{\surd N}\exp({\rm
i}n\phi_{m})=\frac{1}{\surd N}\exp\bigg{(}\frac{2\pi{\rm
i}mn}{N}\bigg{)}=\Phi^{{\scriptscriptstyle(N)}}_{n\pm N}(\phi_{m}),$ (10)

labeled by the cyclic subset $n\in\\{0,1,\ldots,N{-}1\\}$, that are also
orthonormal under the common inner product (8) for discrete functions on
${\cal S}^{1}_{{\scriptscriptstyle(N)}}$, and complete:

$(\Phi^{{\scriptscriptstyle(N)}}_{n},\Phi^{{\scriptscriptstyle(N)}}_{n^{\prime}})_{{\scriptscriptstyle(N)}}=\delta_{n,n^{\prime}},\qquad\sum_{n=0}^{N-1}\\!\Phi^{{\scriptscriptstyle(N)}}_{n}(\phi_{m})^{*}\,\Phi^{{\scriptscriptstyle(N)}}_{n}(\phi_{m^{\prime}})=\delta_{m,m^{\prime}}.$
(11)

Returning to (5) with $(x,y)$ in the polar coordinates $(r,\theta)$ of (3),
and taking for ${f}_{\circ}(\phi_{m})$ the basis functions (10) on the
discrete points of ${\cal S}^{1}_{{\scriptscriptstyle(N)}}$, we write the $N$
solutions to the discretized Helmholtz equation, labeled by cyclical
$n\in\\{0,1,\ldots,N{-}1\\}$, as

$\begin{array}[]{rcl}{}\\!\\!\\!\\!\\!f_{n}(r,\theta_{k})\\!\\!\\!&=&\\!\\!\\!\displaystyle\frac{1}{\surd
N}\sum_{m\in{\cal S}^{1}_{{\scriptscriptstyle(N)}}}\\!\\!\\!\exp[{\rm i}\kappa
r(\cos\theta_{k}\cos\phi_{m}+\sin\theta_{k}\sin\phi_{m})]\,\Phi^{{\scriptscriptstyle(N)}}_{n}(\phi_{m})\\\
{}\\!\\!\\!&=&\\!\\!\\!\displaystyle\frac{1}{N}\sum_{m\in{\cal
S}^{1}_{(N)}}\\!\\!\\!\exp[{\rm i}\kappa r\cos(\theta_{k}-\phi_{m})]\exp({\rm
i}n\phi_{m})\\\ {}\\!\\!\\!&=&\\!\\!\\!\displaystyle\frac{e^{{\rm
i}n(\theta_{k}+\pi/2)}}{N}\sum_{m\in{\cal
S}^{1}_{{\scriptscriptstyle(N)}}}\\!\\!\\!\exp({\rm i}\kappa
r\sin\varphi_{m})\exp(-{\rm i}n\varphi_{m}),\end{array}$ (12)

having replaced $\varphi_{m}:=\theta_{k}-\phi_{m}+\frac{1}{2}\pi$ in the
summation over the $N$ discrete points on the circle.

Following Miller [1, p. 29], the phase in front of (12), $e^{{\rm
i}n\theta_{n}(\theta_{k}+\pi/2)}={\rm i}^{n}e^{2\pi{\rm i}nk/N}={\rm
i}^{n}\sqrt{N}\Phi^{{\scriptscriptstyle(N)}}_{n}(\theta_{k})$, is extracted to
write the functions as

$f_{n}(r,\theta_{k})={\rm
i}^{n}\sqrt{N}\,B^{{\scriptscriptstyle(N)}}_{n}(\kappa
r)\,\Phi^{{\scriptscriptstyle(N)}}_{n}(\theta_{k}),$ (13)

where the radial factor $B^{{\scriptscriptstyle(N)}}_{n}(\rho)$, $\rho:=\kappa
r$, are the discrete Bessel functions. From (12) these functions are seen to
be real and their parities, using coefficients $\\{c_{n},s_{n}\\}:=\\{1,0\\}$
for $n$ even or $\\{0,1\\}$ for $n$ odd, can be written as

$\begin{array}[]{rcl}B^{{\scriptscriptstyle(N)}}_{n}(\rho)&=&\displaystyle\frac{1}{N}\sum_{m\in{\cal
S}^{1}_{N}}\\!\\!\\!\exp({\rm
i}\rho\sin\varphi_{m})\,[c_{n}\cos(n\varphi_{m})-{\rm
i}s_{n}\sin(n\varphi_{m})]\\\ &=&\displaystyle\frac{1}{N}\sum_{m\in{\cal
S}^{1}_{N}}\\!\\!\\!\exp({\rm
i}\rho\sin\varphi_{m})\times\left\\{\begin{array}[]{rl}\cos
n\varphi_{m},&n\hbox{ even},\\\ -{\rm i}\sin n\varphi_{m},&n\hbox{
odd}.\end{array}\right.\end{array}$ (14)

The distinction between even and odd cases of $n$, as done in [3], is subtle
but important to obtain the correct result for all $n$’s (cf. [2, Eq. (9)]).
It results in the parity and cyclicity properties

$B^{{\scriptscriptstyle(N)}}_{n}(\rho)=B^{{\scriptscriptstyle(N)}}_{n\pm
N}(\rho)=(-1)^{n}B^{{\scriptscriptstyle(N)}}_{-n}(\rho)=(-1)^{n}B^{{\scriptscriptstyle(N)}}_{n}(-\rho),\quad
B^{{\scriptscriptstyle(N)}}_{n}(0)=\delta_{n,0},$ (15)

which also hold for the continuous Bessel functions $J_{n}(\rho)$ of integer
order [8].

A plane wave of wavenumber $\kappa$ along the $y$-axis in a Helmholtz medium
that allows only $N$ equidistant directions of propagation on the circle, can
be obtained from (14) using the completeness relation (11) to expand the
middle term and write

$\begin{array}[]{rcl}\exp({\rm
i}\rho\sin\varphi_{m})=B^{{\scriptscriptstyle(N)}}_{0}(\rho)\\!\\!\\!&+&\\!\\!\\!\displaystyle
2\sum_{n=1}^{N-1}B^{{\scriptscriptstyle(N)}}_{2n}(\rho)\cos(2n\varphi_{m})\\\
{}\\!\\!\\!&+&\\!\\!\\!\displaystyle 2{\rm
i}\sum_{n=0}^{N-1}B^{{\scriptscriptstyle(N)}}_{2n+1}(\rho)\sin((2n{+}1)\varphi_{m}),\end{array}$
(16)

showing how the discrete Bessel functions can take the place of the continuous
ones, cf. [7, Eq. KU120(13)].

In Ref. [3] we proved analytically, and verified numerically, that the
following expressions for the discrete Bessel functions are exact analogues of
those valid for continuous Bessel functions. Corresponding to [7, WA44] for
odd $N=:2j+1$, in Ref. [3] we proved the linear relations involving the even
and odd-$n$ discrete Bessel functions,

$\displaystyle B_{0}(\rho)+\sum_{n=1}^{j}B_{2n}(\rho)\cos(2n\varphi_{m})$
$\displaystyle=$ $\displaystyle\cos(\rho\sin\varphi_{m}),$ (17)
$\displaystyle\sum_{n=0}^{j}B_{2n+1}(\rho)\sin((2n{+}1)\varphi_{m})$
$\displaystyle=$
$\displaystyle{\textstyle\frac{1}{2}}\sin(\rho\sin\varphi_{m}).$ (18)

The quadratic formulas [9, §7.6.2, Eq. (6)] associated to the name of Graf,
were shown in Ref. [10] to derive from the rotation of spherical harmonics
through Wigner-$D$ functions, under contraction from the rotation to the
Euclidean group. These relations, of group-theoretical origin, retain their
validity under the discretization of the rotation subgroup, and lead to

$\sum_{n=-2j}^{2j}B_{n}(\rho)\,B_{n^{\prime}-n}(\rho^{\prime})=B_{n^{\prime}}(\rho+\rho^{\prime}),$
(19)

keeping in mind the parity property (15) for the negative $n$-indices in the
sum for odd $N$, addressing the vector rather than spin representations of the
rotation group.

In Fig. 1 we essentially repeat the figure in Ref. [3] where we compared the
discrete and continuous Bessel functions,
$B^{{\scriptscriptstyle(N)}}_{n}(\rho)$ and $J_{n}(\rho)$, to support the
claim that the approximation is indeed remarkable within an interval that is
roughly $0\leq n+\rho<N$. A similar set of figures is presented below for
Mathiew functions.

Figure 1: The ‘discrete’ Bessel functions
$B_{n}^{{\scriptscriptstyle(N)}}(\rho)$ on continuous intervals
$0\leq\rho\leq(2N{-}1)$ (gray lines), vs. the ‘continuous’ Bessel functions
$J_{n}(\rho)$ (thin black lines), for orders $n\in\\{0,\,10,\,30,\,50\\}$ and
point numbers $N\in\\{21,\,61,\,101\\}$. Heavy black lines replace both where
the ‘discrete’ and the ‘continuous’ Bessel functions differ by less than
$10^{-16}$.

Now, having $N$ basis functions $B^{{\scriptscriptstyle(N)}}_{n}(\rho)$,
numbered by cyclic $n$ modulo $N$, it is natural to inquire whether the
argument $\rho$ can or should be also discretized to the $N$ integer values
$\rho_{k}=k\in\\{0,1,\ldots,N{-}1\\}$. This was done in Ref. [2] while in [3]
the plot in Fig. 1 marked these points and used them to define a kernel
$B^{{\scriptscriptstyle(N)}}_{n}(\rho_{k})$ for a ‘discrete Bessel transform’
between two $N$-vectors of components $f_{n}$ and $\widetilde{f}_{k}$. The
fact is that while the angle $\varphi$ is discretized naturally to $N$ points
on the circle, the radial coordinate $\rho$ is not subject to a similarly
compelling set of points, but is valid and non-cyclic over the complex
$\rho$-plane. The same discretization process for the angular —but not the
radial— coordinate will be applied to the Mathieu case below.

##  3 Discrete Mathieu functions

Elliptic coordinates on the plane generalize the previous polar coordinates
(3); they are defined in terms of Cartesian coordinates through

$x=\cosh\varrho\,\cos\psi,\quad
y=\sinh\varrho\,\sin\psi,\qquad\varrho\in[0,\infty),\quad\psi\in(-\pi,\pi]={\cal
S}^{1}.$ (20)

where $(\varrho,\psi)$ are analogues of the previous polar coordinates
$(r,\phi)$ for which we retain the names as ‘radial’ and ‘angular’ variables.
For fixed $\varrho$ or for fixed $\psi$, the locus of points $(x,y)\in{\cal
R}^{2}$ that satisfy

${x^{2}}/{\cosh^{2}\\!\varrho}+{y^{2}}/{\sinh^{2}\\!\varrho}=1,\qquad{x^{2}}/{\cos^{2}\\!\psi}-{y^{2}}/{\sin^{2}\\!\psi}=1,$
(21)

draw families of confocal ellipses or hyperbolas respectively. At $\varrho=0$,
$\psi\in{\cal S}^{1}$ draws twice the line between the two foci $(x,y)=(\pm
1,0)$ for $\psi=(0,\pi)$. The major and minor semi-axes of the ellipses are
$\cosh\varrho$ and $\sinh\varrho$ respectively, so their eccentricities are
$1/\cosh\varrho$, that tend to circles when $\varrho\to\infty$. On the other
hand, for fixed $\psi\in{\cal S}^{1}$ in each of the four quadrants, since
$\varrho\geq 0$, only one of the four arms of the hyperbola is traversed. Thus
we expect four parity cases out of the two reflections, across the $x$ and $y$
axes. Compare this with the case of polar coordinates where $r\geq 0$ but all
reflection axes are equivalent, so $(-1)^{n}$ in (14) provides the two Bessel
parity cases.

Figure 2: Set of equally-spaced discrete points on ellipses (20) of the
‘angular’ coordinates $\\{\psi_{m}\\}\in{\cal
S}^{1}_{{\scriptscriptstyle(N)}}$ for $N=21$, and hyperbolas of the ‘radial’
coordinate for $\varrho\in\\{0.5,\,1,\,1.5\\}$.

The Helmholtz differential equation (2), written in the elliptic coordinates
(20), is clearly separable,

$[(\partial^{2}_{\varrho}+\kappa^{2}\cosh^{2}\\!\varrho)+(\partial^{2}_{\psi}-\kappa^{2}\cos^{2}\\!\psi)]\,f(\varrho,\psi)=0,$
(22)

so that solutions can be written in the product form $f(\varrho,\psi)\sim
P(\varrho)\,\Psi(\psi)$. Dividing by $f$, one obtains two coupled equations in
$\varrho$ and $\psi$, the latter is an eigenvalue equation in the angular
coordinate,

$(\partial^{2}_{\psi}-2q\cos 2\psi)\,\Psi(\psi,q)=\nu\,\Psi(\psi,q),\qquad
q:={\textstyle\frac{1}{4}}\kappa^{2},$ (23)

known as the Mathieu differential equation. The angular coordinate $\psi$ is
periodic and a well-known solution method consists in expanding solutions of
(23) in the Fourier basis $\sim\exp({\rm i}n\psi)$ over all integer $n$. This
defines the Mathieu functions of the first kind ${\rm ce}_{n}(\psi,q)$ and
${\rm se}_{n}(\psi,q)$ with integer $n$ [11], characterized by a parity index
$p\in\\{0,1\\}$ for even and odd cases [7, Eqs. 8.61], and distinct for even
and odd indices. In a two-line expression all cases can be written as

$\bigg{[}{{\rm ce}_{2n+p}(\psi,q)\atop{\rm
se}_{2n+p}(\psi,q)}\bigg{]}=\sum_{s=0}^{\infty}\bigg{[}{A^{2n+p}_{2s+p}\cos((2s{+}p)\psi)\atop
B^{2n+p}_{2s+p}\sin((2s{+}p)\psi)}\bigg{]}.$ (24)

The parities are even ${\rm ce}_{n}(-\psi,q)={\rm ce}_{n}(\psi,q)$, odd ${\rm
se}_{n}(-\psi,q)=-{\rm se}_{n}(\psi,q)$, and ${\rm se}_{0}(\psi,q)\equiv 0$.
The coefficients $A^{n}_{s},\,B^{n}_{s}$ are found introducing this expansion
into (23) to find recursion relations [7, Eqs. 8.62] that lead to efficient
numerical computation. For use below, we write them using Fourier series as

$\bigg{[}{A^{n}_{s}\atop B^{n}_{s}}\bigg{]}=\frac{1}{\pi}\int_{{\cal
S}^{1}}{\rm d}\psi\,\bigg{[}{\cos(s\psi)\,{\rm
ce}_{n}(\psi,q)\atop\sin(s\psi)\,{\rm se}_{n}(\psi,q)}\bigg{]},$ (25)

for $n\neq 0$, while $A^{n}_{0}=(2\pi)^{-1}\int_{{\scriptscriptstyle{\cal
S}^{1}}}{\rm d}\psi\,{\rm ce}_{n}(\psi,q)$, and $B^{n}_{0}\equiv 0$. The
Mathieu functions (24) are orthogonal under the inner product (6) over the
circle, $({\rm ce}_{m},{\rm ce}_{n})_{\circ}=\pi\delta_{m,n}$, $({\rm
se}_{m},{\rm se}_{n})_{\circ}=\pi\delta_{m,n}$ for $n\neq 0$ —zero otherwise,
and $({\rm ce}_{m},{\rm se}_{n})_{\circ}=0$.

We now restrict the range of the angular coordinate $\psi$ from ${\cal S}^{1}$
to ${\cal S}^{1}_{{\scriptscriptstyle(N)}}$, shown for the elliptic
coordinates in Fig. 2, in correspondence with the previous discrete phase
functions in the Bessel case (10), and thus defining the ‘angular’ discrete
Mathieu functions of the first type over $\psi_{m}\in{\cal
S}^{1}_{{\scriptscriptstyle(N)}}$ as

$\bigg{[}{{\rm ce}^{{\scriptscriptstyle(N)}}_{2n+p}(\psi_{m},q)\atop{\rm
se}^{{\scriptscriptstyle(N)}}_{2n+p}(\psi_{m},q)}\bigg{]}:=\sum_{s=0}^{N-1}\bigg{[}{a^{2n+p}_{2s+p}\cos((2s{+}p)\psi_{m})\atop
b^{2n+p}_{2s+p}\sin((2s{+}p)\psi_{m})}\bigg{]},$ (26)

with coefficients $a^{n}_{s},\,b^{n}_{s}$. The finite $N$-point Fourier
transform approximates them through the replacement (7) to the functions and
coefficients $A^{n}_{s},\,B^{n}_{s}$ of the continuous case in (25), as

$\bigg{[}{a^{n}_{s}\atop
b^{n}_{s}}\bigg{]}:=\frac{1}{N}\sum_{m=0}^{N-1}\bigg{[}{\cos(s\psi_{m})\,{\rm
ce}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)\atop\sin(s\psi_{m})\,{\rm
se}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)}\bigg{]}\simeq\frac{1}{2}\bigg{[}{A^{n}_{s}\atop
B^{n}_{s}}\bigg{]},$ (27)

for $s\neq 0$, while $a_{0}^{n}=A_{0}^{n}$, $b^{n}_{0}=0$, and also
$b^{0}_{s}=0$.

Figure 3: Discrete vs. continuous ‘angular’ Mathieu functions for $N=41$,
$q=2$. The values of the discrete functions ${\rm
ce}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)$ and ${\rm
se}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)$ at $\psi_{m}$, $0\leq m\leq
N{-}1$, are indicated by circles. The continuous Mathieu functions ${\rm
ce}_{n}(\psi,q)$ and ${\rm se}_{n}(\psi,q)$ are marked by black lines in their
full range $0\leq\psi<2\pi$. Their difference is less than $10^{-16}$ for all
points $\psi_{m}$.

The last relation in (27) is an approximate equality, the validity of which is
contingent upon the numerical computation and comparison between the lower-
and upper-case coefficients within a range of their indices in, say, $0\leq
n,s\leq N{-}1$, which is reflected in turn by the discrete and continuous
Mathieu functions themselves. In Fig. 3 we compare a sample of continuous
angular Mathieu functions of the first kind with their discrete approximations
from Eq. (26). In favor of the thus defined discrete Mathieu functions, we
note that they satisfy orthogonality relations under the discrete inner
product (8), namely

$({\rm ce}^{{\scriptscriptstyle(N)}}_{n},{\rm
ce}^{{\scriptscriptstyle(N)}}_{n^{\prime}})_{{\scriptscriptstyle(N)}}={\textstyle\frac{1}{2}}N\delta_{n,{n^{\prime}}},\\!\quad({\rm
se}^{{\scriptscriptstyle(N)}}_{n},{\rm
se}^{{\scriptscriptstyle(N)}}_{{n^{\prime}}\neq
0})_{{\scriptscriptstyle(N)}}={\textstyle\frac{1}{2}}N\delta_{n,{n^{\prime}}},\\!\quad({\rm
ce}^{{\scriptscriptstyle(N)}}_{n},{\rm
se}^{{\scriptscriptstyle(N)}}_{n^{\prime}})_{{\scriptscriptstyle(N)}}=0.$ (28)

By construction, the parities of the discrete Mathieu functions are also even
${\rm ce}^{{\scriptscriptstyle(N)}}_{n}(-\psi_{m},q)={\rm
ce}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)$, or odd ${\rm
se}^{{\scriptscriptstyle(N)}}_{n}(-\psi_{m},q)=-{\rm
se}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)$.

At this point it is illuminating to inquire into the manner in which the
discrete functions approximate the continuous ones. Consider for example how
${\rm ce}^{{\scriptscriptstyle(N)}}_{0}(\psi_{m},q)$, whose definition (26)
allows us to compute it for continuous $\psi_{m}\in{\cal S}^{1}$, matches
${\rm ce}_{0}(\psi,q)$ in the whole $\psi$ range. In Fig. 4 we do so for small
$N$, noting that where the continued $\psi_{m}$ lines of the former take their
values, they intersect the properly continuous line of the latter; although
the two lines intersect also at other points, the two lines remain notably
distinct. As the figure shows, the approximation is not valid over presumably
small ranges around these intersections, but only at the prescribed
$\psi_{m}=2\pi m/N$ points. We intend to elaborate on such and similar limits
elsewhere.

Figure 4: Comparison between the ‘discrete’ Mathieu functions ${\rm
ce}^{{\scriptscriptstyle(N)}}_{0}(\psi_{m},q)$ whose arguments are continued
to $\psi_{m}\in{\cal S}^{1}$ (gray line) vs. the ‘continuous’ Mathieu function
${\rm ce}_{n}(\psi,q)$ (black line), for point numbers $N\in\\{5,11,21\\}$ and
$q=2$. The discrete points $\psi_{m}\in{\cal S}^{1}_{{\scriptscriptstyle(N)}}$
lie at a subset of the intersections marked with circles.

Proceeding now as we did in (12), but using the discrete Mathieu functions of
the first kind ${\rm ce}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)$ and ${\rm
se}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)$ in place of the plain phase
functions $\Phi^{{\scriptscriptstyle(N)}}_{n}(\phi_{m})$, we again have
Helmholtz solutions $f_{n}(\varrho,\psi_{k})$ on the plane that are
characterized by parities and sub-indices $n$, whose radial factor will be the
discrete ‘radial’ Mathieu functions of the second kind, to be indicated
correspondingly as ${\rm Ce}^{{\scriptscriptstyle(N)}}_{n}(\varrho,q)$ and
${\rm Se}^{{\scriptscriptstyle(N)}}_{n}(\varrho,q)$,

$\displaystyle{}\\!\\!\\!\\!\\!\\!\\!\\!\\!\bigg{[}{f^{{\scriptscriptstyle
c}}_{2n+p}(\varrho,\psi_{k})\atop f^{{\scriptscriptstyle
s}}_{2n+p+1}(\varrho,\psi_{k})}\bigg{]}$ $\displaystyle=$
$\displaystyle\frac{1}{N}\sum_{m=1}^{N-1}\bigg{[}{{\rm
ce}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)\atop{\rm
se}^{{\scriptscriptstyle(N)}}_{n}(\psi_{m},q)}\bigg{]}\exp[{\rm
i}\kappa(x\cos\psi_{m}+y\sin\psi_{m})$ (29) $\displaystyle=:$
$\displaystyle\bigg{[}{c_{n}(q)\,{\rm
Ce}^{{\scriptscriptstyle(N)}}_{n}(\varrho,q)\,{\rm
ce}^{{\scriptscriptstyle(N)}}_{n}(\psi_{k},q)\atop s_{n}(q)\,{\rm
Se}^{{\scriptscriptstyle(N)}}_{n}(\varrho,q)\,{\rm
se}^{{\scriptscriptstyle(N)}}_{n}(\psi_{k},q)}\bigg{]},$ (30)

where $c_{n}(q)$ and $s_{n}(q)$ are constants. Using the elliptic coordinates
with a discretized angular part,
$x(\varrho,\psi_{k})=\cosh\varrho\\*\cos\psi_{k}$ and
$y(\varrho,\psi_{k})=\sinh\varrho\\*\sin\psi_{k}$ in (20), the phase exponent
is then ${\rm i}\kappa=2{\rm i}\surd q$ times
$x\cos\psi_{m}+y\sin\psi_{m}=\cosh\varrho\\*\cos\psi_{k}\\*\cos\psi_{m}+\sinh\varrho\\*\sin\psi_{k}\\*\sin\psi_{m}$.
As was done before in (12) and (13), we extract the new discrete ‘radial’
functions using the orthogonality (28) of the previous discrete ‘angular’
Mathieu functions, as

${}\\!\\!\\!\\!\\!\begin{array}[]{r}\displaystyle\bigg{[}{{\rm
Ce}^{{\scriptscriptstyle(N)}}_{2n+p}(\varrho,q)\atop{\rm
Se}^{{\scriptscriptstyle(N)}}_{2n+p+1}(\varrho,q)}\bigg{]}=\displaystyle\bigg{[}{1/N\,{c}_{2n+p}(q)\,{\rm
ce}^{{\scriptscriptstyle(N)}}_{2n+p}(\psi_{k},q)\atop
1/N\,{s}_{2n+p+1}(q)\,{\rm
se}^{{\scriptscriptstyle(N)}}_{2n+p+1}(\psi_{k},q)}\bigg{]}\sum_{m=0}^{N-1}\bigg{[}{{\rm
ce}^{{\scriptscriptstyle(N)}}_{2n+p}(\psi_{m},q)\atop{\rm
se}^{{\scriptscriptstyle(N)}}_{2n+p+1}(\psi_{m},q)}\bigg{]}\\\\[10.0pt]
{}\times\exp[2{\rm
i}\sqrt{q}\,(\cosh\varrho\\*\cos\psi_{k}\\*\cos\psi_{m}+\sinh\varrho\\*\sin\psi_{k}\\*\sin\psi_{m})].\end{array}$
(31)

The coefficients in front of the summation will be now determined through
considering specific values for the ‘angular’ coordinate
$\psi\leftrightarrow\psi_{m}$, comparing them with expressions of the
continuous Mathieu functions of the second kind obtained from integrals that
are tabulated in Ref. [7, §6.92]. There, the exponential factors appear with
only a single summand in the exponent, either sine or cosine. This occurs in
(31) only for $\psi_{m}=0$ or ${\textstyle\frac{1}{2}}\pi$, although the
latter is not in the set ${\cal S}_{{\scriptscriptstyle(N)}}^{1}$ if
$\psi_{0}=0$, since $N$ was assumed to be odd.

Let us first consider the case of even parity $p=0$ and the angle
$\psi={\textstyle\frac{1}{2}}\pi$ in (31), where the previous remark applies.
Based on the close approximation between the discrete and continuous Mathieu
functions, we may simply replace the latter for the former, so that the two
lines in that expression read

$\bigg{[}{{\rm Ce}^{{\scriptscriptstyle(N)}}_{2n}(\varrho,q)\atop{\rm
Se}^{{\scriptscriptstyle(N)}}_{2n+1}(\varrho,q)}\bigg{]}=\bigg{[}{K^{{\scriptscriptstyle
c}}_{2n}\atop K^{{\scriptscriptstyle
s}}_{2n+1}}\bigg{]}\sum_{m=0}^{N-1}\bigg{[}{{\rm
ce}^{{\scriptscriptstyle(N)}}_{2n}(\psi_{m},q)\atop{\rm
se}^{{\scriptscriptstyle(N)}}_{2n+1}(\psi_{m},q)}\bigg{]}\exp(2{\rm
i}\sqrt{q}\,\sinh\varrho\sin\psi_{m}).$ (32)

When this summation formula is compared with the integral expressions
tabulated in [7, §6.92], namely

$\begin{array}[]{r}\displaystyle\bigg{[}{{\rm Ce}_{2n}(\varrho,q)\atop{\rm
Se}_{2n+1}(\varrho,q)}\bigg{]}=\bigg{[}{{\rm
ce}_{2n}(0,q)/2\pi\,A_{0}^{2n}\atop-{\rm i}{\rm
se}^{\prime}_{2n+1}(0,q)/2\pi\,B_{1}^{2n+1}\\!\\!\sqrt{q}}\bigg{]}{}\\\\[10.0pt]
\displaystyle{}\times\int_{{\cal S}^{1}}{\rm d}\psi\,\bigg{[}{{\rm
ce}_{2n}(\psi,q)\atop{\rm se}_{2n+1}(\psi,q)}\bigg{]}\exp(2{\rm
i}\sqrt{q}\,\sinh\varrho\sin\psi),\end{array}$ (33)

we conclude that the constants in the summation (32), after identifying
$2\pi\leftrightarrow N$, $A_{0}^{2n}=a_{0}^{2n}$ and $B_{1}^{2n+1}\simeq
2b_{1}^{2n+1}$, are

$K^{{\scriptscriptstyle c}}_{2n}=\frac{{\rm ce}_{2n}(0,q)}{a_{0}^{2n}N},\qquad
K^{{\scriptscriptstyle s}}_{2n+1}=\frac{-{\rm i}\,{\rm
se}^{\prime}_{2n+1}(0,q)}{2b_{1}^{2n+1}N\surd q}.$ (34)

where ${\rm se}^{\prime}_{n}(0,q):={\rm d}\,{\rm se}_{n}(\psi,q)/{\rm
d}\psi|_{\psi=0}$. In Fig. 5 we compare a sample of the discrete and
continuous ‘radial’ Mathieu functions, noting that the two lines are quite
coincident in the range $\varrho\in[0,\pi)$, but that the discrete approximant
oscillates wildly beyond $\pi$. Again, here we can only justify this statement
numerically.

Next we consider the case of odd parity $p=1$ at the value $\psi=\psi_{0}=0$.
The upper line in (32) reads

${\rm Ce}^{{\scriptscriptstyle(N)}}_{2n+1}(\varrho,q)=K^{{\scriptscriptstyle
c}}_{2n+1}\sum_{m=0}^{N-1}{\rm
ce}^{{\scriptscriptstyle(N)}}_{2n+1}(\psi_{m},q)\exp(2{\rm
i}\sqrt{q}\cosh\varrho\cos\psi_{m}),$ (35)

that we compare with the integral for the continuous Mathieu functions of the
second kind in [7, §6.92], namely

${\rm Ce}_{2n+1}(\varrho,q)=\frac{{\rm i}\,{\rm
ce}^{\prime}_{2n+1}({\textstyle\frac{1}{2}}\pi,q)}{2\pi
A_{1}^{2n+1}\sqrt{q}}\int_{{\cal S}^{1}}{\rm d}\psi\,{\rm
ce}_{2n+1}(\psi,q)\exp(2{\rm i}\sqrt{q}\cosh\varrho\cos\psi),$ (36)

where ${\rm ce}^{\prime}_{2n+1}(\psi,q)$ is the derivative of the Mathieu
function. Again exploiting the correspondences (7), $2\pi\leftrightarrow N$
and $A_{1}^{2n+1}\simeq 2a_{1}^{2n+1}$, we conclude the constant in (38) to be

$K^{{\scriptscriptstyle c}}_{2n+1}=\frac{{\rm i}\,{\rm
ce}^{\prime}_{2n+1}({\textstyle\frac{1}{2}}\pi,q)}{2a_{1}^{2n+1}N\surd q}.$
(37)

Figure 5: Discrete vs. continuous ‘radial’ Mathieu functions in the interval
$0\leq\varrho<3.3$, for $N\in\\{5,11,21\\}$ and here for $q=2$. The ‘discrete’
functions ${\rm Ce}^{{\scriptscriptstyle(N)}}_{n}(\varrho,q)$ and ${\rm
Se}^{{\scriptscriptstyle(N)}}_{n}(\varrho,q)$ with the (continuous) argument
$\varrho$ (gray line), is compared with the ‘continuous’ functions ${\rm
Ce}_{n}(\varrho,q)$ and ${\rm Se}_{n}(\varrho,q)$ (thin black line). As
before, where both coincide within $10^{-16}$ they are replaced by a thick
black line. The radial Mathieu functions, when computed with the commercial
Mathematica algorithm, oscillate wildly after an upper value that decreases
with increasing values of $q$.

The remaining case to be considered is that of odd parity and even index,
namely for ${\rm Se}^{{\scriptscriptstyle(N)}}_{2n+2}(\varrho,q)$. This
presents a problem though, because the summation (31) is identically zero for
both $\psi_{k}=0$ and ${\textstyle\frac{1}{2}}\pi$ due to the parities of the
terms in the sum. It is different from zero for
$0<\psi_{k}<{\textstyle\frac{1}{2}}\pi$ however, so if we choose
$\psi_{k}=\frac{1}{4}\pi$, where both summands in the exponent appear as
$1/\surd 2$, we can write

$\begin{array}[]{r}\displaystyle{\rm
Se}^{{\scriptscriptstyle(N)}}_{2n+2}(\varrho,q)=K^{{\scriptscriptstyle
s}}_{2n+2}\sum_{m=0}^{N-1}{\rm
se}^{{\scriptscriptstyle(N)}}_{2n+2}(\psi_{m},q){}\\\\[5.0pt]
\displaystyle{}\times\exp[{\rm
i}\sqrt{2q}\,(\cosh\varrho\cos\psi_{m}+\sinh\varrho\sin\psi_{m})].\end{array}$
(38)

For the corresponding continuous case, we could not find a corresponding
integral in [7, §6.92], so we cannot give a closed expression for the
coefficient $K^{{\scriptscriptstyle s}}_{2n+2}$ in (38). The lack of a similar
plane-wave integral expression for the continuous Mathieu functions ${\rm
Se}_{2n+2}(\varrho,q)$ has been noted also in Ref. [12] without explanation.
However, we have checked numerically that the simile of the discrete to
continuous functions approximation provided by

${\rm Se}^{{\scriptscriptstyle(N)}}_{2n+2}(\varrho,q)\simeq-{\rm i}\,{\rm
se}_{2n+2}({\rm i}\varrho,q)={\rm Se}_{2n+2}(\varrho,q),$ (39)

which is an equality for continuous functions, cf. [7, 8.611.4, 8.631.4]. For
$0<\varrho<2$ the difference in (39) less than $10^{-14}$. We should note that
generally the discrete ‘radial’ and ‘angular’ Mathieu functions for pure
imaginary arguments are not related to similar equalities of their continuous
integral expressions, because the summation definitions in (26) involve
hyperbolic functions. In particular, say,

${\rm Se}^{{\scriptscriptstyle(N)}}_{2n+p+1}(\varrho,q)\neq-{\rm i}{\rm
se}^{{\scriptscriptstyle(N)}}_{2n+p+1}({\rm
i}\varrho,q)=\sum_{m=0}^{N-1}b^{2n+p+1}_{2m+p+1}\sinh[(2m{+}p{+}1)\varrho].$
(40)

##  4 Concluding remarks

The expansion of Helmholtz plane waves in series of radial Bessel and angular
trigonometric functions has its discrete analogue in Eq. (16), which tells us
that the wavefield due to a finite number $N$ of plane waves at equidistant
direction angles can be expanded in discrete Bessel radial functions and
corresponding trigonometric angular functions. A similar statement will hold
when the wavefield is expanded in discrete Mathieu functions with the phases
determined by the points on an ellipse as depicted in Fig. 2. Conceivably such
fields can be produced in resonant two-dimensional micro-cavities fed by a
number of activation channels.

We recognize that the full treatment and exploration of properties for the
discrete Bessel and Mathieu function presented here is not exhaustive, but
that it should be sufficient to indicate that the approximation method
consisting in the replacement of a continuous closed subgroup of the symmetry
group of a partial differential equation by a finite discrete group is
definitely of interest. In the present case of two dimensions, the orthogonal
group was reduced to the dihedral group. In three dimensions, the symmetry
Euclidean symmetry group could reduce its three-dimensional rotation subgroup
by any of its polyhedral subgroups, whose functions may serve to describe
wavefields with a corresponding subset of wave propagation directions. Here we
have presented a set of exact relations, others whose approximation closeness
was estimated through numerical computation, and others that have been only
suggested by that approach, and for which we expect to present further results
from ongoing work.

## Acknowledgments

We thank the support of the Universidad Nacional Autónoma de México through
the PAPIIT-DGAPA project AG100120 Óptica Matemática.

## References

  * [1] W. Miller Jr., Symmetry and Separation of Variables, Encyclopedia of Mathematics, Vol. 4 (Cambridge Univerity Press, 1984).
  * [2] G. Biagetti, P. Crippa, L. Falaschetti, and C. Turchetti, Discrete Bessel functions for representing the Class of Finite Duration Decaying Sequences, European Signal Analysis Conference, pp. 2126–2130 (Budapest, 2016).
  * [3] K. Uriostegui and K.B. Wolf, Discrete Bessel functions and transform,
(submitted) arXiv:2005.06076 [math-ph]

  * [4] R.H. Boyer, Discrete Bessel functions, J. Math. Anal. Appl. 2, 509–524 (1961).
  * [5] M. Bohner and T. Cuchta, The Bessel difference equation, Proc. Amer. Math. Soc. 145, 1567–1580 (2017).
  * [6] A. Slavík, Discrete Bessel functions and partial differential equations, J. Diff. Eqs. Applics. DOI:10.1080/10236198.2017.141610 (2017).
  * [7] I.S. Gradshteyn and I.M. Ryzhik, Tables of Integrals, Series and Products, A. Jeffrey and D. Zwillinger Eds, (Academic Press, 2007).
  * [8] G.N. Watson, Theory of Bessel Functions (Cambridge University Press, 1922).
  * [9] A. Erdélyi et al., Higher Transcendental Functions (Based on notes by H. Bateman) Vol. 2 (McGraw-Hill, New York, 1953).
  * [10] P. Winternitz, K.B. Wolf, G.S. Pogosyan, and A.N. Sissakian, Graf’s addition theorem obtained from SO(3) contraction, Theor. Mat. Phys. 129, 1501–1503 (2001).
  * [11] N.W. McLachlan, Theory and Application of Mathieu Functions (Oxford University Press, 1947).
  * [12] L. Chaos-Cador and E. Ley-Koo, Mathieu functions, matrix evaluation, and generating functions, Rev. Mex. Fís. 48, 67–75 (2002).


and perhaps these effects are correlated with the bulk properties of the three
CV classes.

B. A surprising open question is to understand the existence of the coherent
non-orbital periods (see Table 4), with these novae each having 1–3 such
periods, with most being perplexing and mysterious in origin. Five of the
novae (V842 Cen, QY Mus, V598 Pup, XX Tau, and PW Vul) have coherent, long-
lasting, and highly significant periodicities from 0.12–4.1 days, which are
certainly not the orbital period. Two of the systems (V2574 Oph and V407 Lup)
have fast periodicities that are stable and coherent, with the only plausible
good ‘clock’ being the white dwarf spin, yet for which they are not the spin
period, and there must be some unknown clocking mechanism. The unique new
phenomena for YZ Ret can be speculated to be associated with some exotic
superhump condition, with no precedent or prior theory. These modulations are
mysterious because their coherence forces them to each be tied to some
accurate ‘clock’ in the system, but they cannot be either orbital or spin,
while no other highly-coherent clocks can be invoked. Throughout astrophysics,
whenever a periodicity is discovered, physics can be applied to the system,
while I would be hopeful that an understanding of these coherent non-orbital
periods can provide the physics of some previously unknown mechanisms
operating in nova binaries. To the best of my knowledge of the literature, no
one has ever recognized or highlighted this basic problem.

C. Novae with evolved companions are now seen to make up 31 per cent of all
Milky Way novae. These novae have been completely ignored by modelers. No idea
has been published as to the past and future changes in the accretion rate,
$\dot{P}$, and stellar masses. Is the accretion being driven by a decrease in
$P$, or by the evolutionary expansion of the companion? Indeed, the dominant
mechanism for angular momentum loss in the binary is not known. A conundrum
recognized in this paper is that the novae with red giant companions are 75
per cent in the bulge population, with this being difficult to understand in
detail, especially when the novae with red giant companions in the Andromeda
Galaxy comprise $\sim$30 per cent of the observed novae and these all appear
to be in the disc population (Williams et al. 2016). Another unaddressed
question is whether the currently observed systems have only recently come
into contact. Why are these novae so sharply different from symbiotic novae?
With some evolved companions likely being more massive than the white dwarf,
will the simple mass transfer drive runaway accretion? Is the white dwarf
gaining or losing mass over each eruption cycle, with this leading into the
question of whether CVs are the solution to the highly-important Type Ia
supernova Progenitor Problem?

##  10 Acknowledgements

I thank Rebekah Hounsell for helping me to get working the Lightkurve program
available from MAST. The AAVSO provides the unique tools that formed the
backbone of my work, including variable star cataloguing (VSX), their huge
International Database of light curves, their finder chart facility (VSP),
their calibration of comparison stars (APASS), their DFT routine (VSTAR), plus
their ability to respond to get new special-purpose data. My use of the high
level data products of the many sky surveys does not well highlight the huge
efforts over many years by large teams of workers from $TESS$, ZTF, K2, MAST,
AAVSO, SMARTS, ASAS, and OGLE. Funding for the $TESS$ mission is provided by
NASA’s Science Mission directorate. This research made use of Lightkurve, a
Python package for $Kepler$ and $TESS$ data analysis (Lightkurve
Collaboration, 2018).

##  11 Data Availability

All of the photometry data used in this paper are publicly available from the
references and links in Section 2.1.

## References

  * Ak A1 (2005) Ak T., Retter A., Liu A., 2005, Publ. Astron. Soc. Australia, 22, 298
  * Aydi et al. (2018) Aydi E. et al., 2018, MNRAS, 480, 572
  * Aydi A2 (2019) Aydi E. et al., 2019, ATel, 12889
  * Aydi A3 (2021) Aydi E. et al., 2021, arXiv:2108.07868
  * Baklanov A4 (2008) Baklanov A., Pavlenko E., Berezina E., 2008, ATel, 1514
  * Balman A5 (2009) Balman S., Nasiroglu I., Akyuz A., 2009, ATel, 2137
  * Balman et al. (2006) Balman S., Retter A., Bos M., 2006, AJ, 131, 2628
  * Balman A6 (2005) Balman S., Yilmaz A., Retter A., Saygac T., Esenoglu H., 2005, MNRAS, 356, 773
  * Baptista A7 (1993) Baptista R., Jablonski F. J., Cieslinski D., Steiner J. E., 1993, ApJ, 406, L67
  * Beardmore A8 (2012) Beardmore A. P. et al., 2012, A&A, 545, A116
  * Bellm et al. (2019) Bellm E. C. et al., 2019, PASP, 131, 18002
  * Bianchini A9 (2012) Bianchini A., Saygac T., Orio M., Della Valle M., Williams R., 2012, A&A, 539, A94
  * Bos et al. (2001a) Bos M., Retter A., McCormick J., Velthuis F., 2001a, IAUCirc, 7610
  * Bos A10 (2001b) Bos M., Retter A., McCormick J., Velthius F., 2001b, IAUCirc, 7665
  * Brandi A11 (2009) Brandi E., Quiroga C., Mikolajewski J., Ferrer O., Garcia L. G., 2009, A&A, 497, 815
  * Chochol A12 (2013) Chochol D., Shugarov S. Y., Volkov I. M., Goranskij V. P., Metlova N. V., Barsukova E. A., Gabdeev M. M., 2013, IBVS, 6045
  * Dai A13 (2010) Dai Z., Qian S., 2010, New Astron., 15, 380
  * Denisenko A14 (2019) Denisenko D., 2019, vsnet-alert, 23463
  * Diaz and Bruch (1994) Diaz M. P., Bruch A., 1994, IBVS, 4079
  * Diaz A15 (1997) Diaz M. P., Bruch A., 1997, A&A, 322, 807
  * Dobrotka A16 (2006a) Dobrotka A., Friedjung M., Retter A., Hric L., Novak R., 2006a, A&A, 448, 1107
  * Dobrotka A17 (2006b) Dobrotka A., Retter A., Liu A., 2006b, MNRAS, 371, 459
  * Downes et al. (2001) Downes R. A., Webbink R. F., Shara M. M., Ritter H., Kolb U., Duerbeck H. W., 2001, PASP, 113, 764
  * Duerbeck (1987) Duerbeck H. W., 1987, Space Science Reviews, 45, 1
  * Dyck et al. (1996) Dyck H. M., Benson J., van Belle G., Ridgway S., 1996, AJ, 111, 1705
  * Egan et al. (2014) Egan J. M., Woudt P. A., Warner B., Williams R. E., Steeghs D. T. H., Ribeiro V. A. R. M., 2014, in Woudt P. A., Ribeiro V. A. R. M., eds, Astronomical Society of the Pacific Conference Series Vol. 490, Stellar Novae: Past and Future Decades. p. 67
  * Fuentes-Morales et al. (2021) Fuentes-Morales I. et al., 2021, MNRAS, 501, 6083
  * Goranskij A18 (2006) Goranskij V. P., Metlova N. V., Burenkov A. N., 2006, ATel, 928
  * Goranskij A19 (2010) Goranskij V., Shugarov S., Zharova A., Kroll P., Barsukova E. A., 2010, Peremennye Zvezdy, 30
  * Goransky A20 (1997) Goransky V. P., Shugarov S. Y., Dmitrienko E. S., Pavlenko E. P., 1997, in Maoz D., Sternberg A., Leibowitz E. M., eds, As- trophysics and Space Science Library Vol. 218, Astronomical Time Series. pp. 219-222
  * Hachisu and Kato (2019) Hachisu I., Kato M., 2019, ApJS, 242, 18
  * Hacke (1987) Hacke G., 1987, Mitt. Verand. Sterne, 11, 40
  * Haefner A21 (1999) Haefner R., 1999, IBVS, 4706
  * Harrison A22 (1992) Harrison T. E., 1992, MNRAS, 259, 17P
  * Hatano et al. (1997) Hatano K., Branch D., Fiser A., Starrfield S., 1997, MNRAS, 290, 113
  * Howell et al. (2014) Howell S. B. et al., 2014, PASP, 126, 398
  * Ingram A23 (1992) Ingram D., Garnavich P., Green P., Szkody P., 1992, PASP, 104, 402
  * Joshi A24 (2015) Joshi V., Banerjee D. P. K., Ashok N. M., Venkataraman V., Walter F. M., 2015, MNRAS, 452, 3696
  * Kaluzny A25 (1990) Kaluzny J., 1990, MNRAS, 245, 547
  * Kang A26 (2006) Kang T. W., Retter A., Liu A., Richards M., 2006, AJ, 132, 608
  * Kato A27 (2004) Kato T., Ishioka R., Uemura M., Starkey D., Krajci T., 2004, PASJ, 56, S125
  * Katysheva A28 (2007) Katysheva N. A., Shugarov S. Y., 2007, in Napiwotzki R., Burleigh M. R., eds, Astronomical Society of the Pacific Con- ference Series Vol. 372, 15th European Workshop on White Dwarfs. p. 523
  * Knigge et al. (2011) Knigge C., Baraffe I., Patterson J., 2011, ApJS, 194, 28
  * Kraft (1964) Kraft R. P., 1964, ApJ, 139, 457
  * Landolt (2016) Landolt, A. U., 2016, JAAVSO, 44, 45
  * Leibowitz A29 (1997) Leibowitz E. M., Ofek E. O., Mattei J. A., 1997, MNRAS, 287, 634
  * Lightkurve Collaboration (2018) Lightkurve Collaboration, 2018, Astrophysics Source Code Library, record ascl:1812.013
  * Marin and Shafter (2009) Marin E., Shafter A. W., 2009, PASP, 121, 1090
  * Mikolajewska (2010) Mikołajewska J., 2010, arXiv:1011.5657
  * Mikolajewska et al. (2021) Mikołajewska J. et al., 2021, MNRAS, 2021, 504, 2122
  * Morales-Rueda A30 (2002) Morales-Rueda L., Still M. D., Roche P., Wood J. H., Lockley J. J., 2002, MNRAS, 329, 597
  * Mroz et al. (2014) Mróz P. et al., 2014, MNRAS, 443, 784
  * Mroz et al. (2015) Mróz P. et al., 2015, ApJS, 219, 26
  * Mroz A32 (2016) Mróz P. et al., 2016, Nature, 537, 649
  * Mroz A31 (2021) Mróz P., Udalski A., OGLE team, 2021, ATel, 14410
  * Munari A33 (2018) Munari U., Banerjee D. P. K., 2018, MNRAS, 475, 508
  * Munari A34 (2013) Munari U., Dallaporta S., Castellani F., Valisa P., Frigo A., Chomiuk L., Ribeiro V. A. R. M., 2013, MNRAS, 435, 771
  * Narloch A35 (2014) Narloch W., Kaluzny J., Krzeminski W., Pych W., Rozyczka M., Shectman S., Thompson I. B., Tomov T., 2014, Baltic Astronomy, 23, 1
  * Ochner A36 (2015) Ochner P., Moschini F., Munari U., Frigo A., 2015, MNRAS, 454, 123
  * Olech A37 (2002) Olech A., 2002, Acta Astron., 52, 273
  * Ozdonmez et al. (2018) Özdönmez A., Ege E., Güver T., Ak T., 2018, MNRAS, 476, 4162
  * Pagnotta and Schaefer (2014) Pagnotta A., Schaefer B. E., 2014, ApJ, 788, 164
  * Pagnotta et al. (2009) Pagnotta A., Schaefer B. Xiao L., Collazzi A., Kroll P., 2009, AJ, 138, 1230
  * Patterson (1984) Patterson J., 1984, ApJS, 54, 443
  * Patterson A41 (2010) Patterson J., Halpern J. P., Stockdale C., Bolt G., Lowther S., Monard B., 2010, ATel, 2746
  * Patterson A39 (2017) Patterson J. O. et al., 2017, MNRAS, 466, 581
  * Patteron A38 (2021b) Patterson J., Epstein-Martin M., Vanmunster T., Kemp J., 2021, ATel, 14856
  * Patterson A40 (2022) Patterson J. O. et al., 2022, ApJ, 924, 27
  * Pavlenko A42 (2018) Pavlenko E. P., et al., 2018, MNRAS, 479, 341
  * Payne-Gaposchkin (1964) Payne-Gaposchkin C. 1964, The Galactic Novae, Dover, New York
  * Peters A43 (2006) Peters C. S., Thorstensen J. R., 2006, PASP, 118, 687
  * Plaut (1958) Plaut L., 1958, Leiden Ann, 21, 3
  * Pojmanski (1997) Pojmanski G., 1997, Act Astron., 47, 467
  * Rappaport, Verbunt, and Joss (1983) Rappaport S., Verbunt F., Joss P. C., 1983, ApJ, 275, 713
  * Retter A44 (1995) Retter A., Leibowitz E. M., 1995, IAUCirc, 6234
  * Retter A45 (1998) Retter A., Leibowitz E. M., Kovo-Kariti O., 1998, MNRAS, 293, 145
  * Retter A46 (1999a) Retter A., Leibowitz E. M., Naylor T., 1999a, MNRAS, 308, 140
  * Retter A47 (1999b) Retter A., Liller W., Garradd G., 1999b, IAU Circ., 7124
  * Ricker et al. (2015) Ricker G. R. et al., 2015, Jour. Astron. Telescopes Instr. Systems, 1, 014003
  * Ringwald A49 (1996) Ringwald F. A., Naylor T., Mukai K., 1996, MNRAS, 281, 192
  * Ringwald A48 (2005) Ringwald F. A., Chase D. W., Reynolds D. S., 2005, PASP, 117, 1223
  * Ritter and Kolb (2003) Ritter H., Kolb U. 2003, A&A, 404, 301 (updated RKcat version 7.24, 2016)
  * Rodriguez-Gil and Torres (2005) Rodríguez-Gil P., Torres M. A. P., 2005, A&A, 431, 289
  * Rodr A50 (2010) Rodríguez-Gil P., et al., 2010, MNRAS, 407, L21
  * Saito et al. (2013) Saito R. K. et al., 2013, A&A, 554, 123
  * Salazar A51 (2017) Salazar I., LeBleu A., Schaefer B., Landolt A., 2017, MNRAS, 469, 4116
  * Samus et al. (2017) Samus N.N., Kazarovets E.V., Durlevich O.V., Kireeva N.N., Pastukhova E.N., 2017, Astronomy Reports, 61, 80
  * Schaefer (2009) Schaefer B. E., 2009, ApJ, 697, 721
  * Schaefer (2010) Schaefer B. E., 2010, ApJS, 187, 275
  * Schaefer A53 (2011) Schaefer B. E., 2011, ApJ, 742, 112
  * Schaefer A54 (2020) Schaefer B. E., 2020, MNRAS, 492, 3323
  * Schaefer (2021b) Schaefer B. E., 2021, RNAAS, 5, 150
  * Schaefer and Collazzi (2010) Schaefer B. E., Collazzi A., 2010, AJ, 139, 1831
  * Schaefer et al. (1992) Schaefer B. E. et al., 1992, ApJS, 81, 321
  * Schaefer A52 (2019) Schaefer B. E. et al., 2019, MNRAS, 487, 1120
  * Schaefer, A57 (2020) Schaefer, B. E., 2020, MNRAS, 492, 3343
  * Schaefer A55 (2022a) Schaefer B. E., Pagnotta A., Zoppelt S., 2022a, MNRAS, 512, 1924
  * Schaefer A56 (2022b) Schaefer B. E., Walter F. M., Hounsell R., Hillman Y., 2022b, MNRAS, submitted
  * Schlafly and Finkbeiner (2011) Schlafly E.F., Finkbeiner D.P., 2011, ApJ 737, 103
  * Schmidt A58 (2020a) Schmidt R. E., 2020a, JAAVSO, 48, 13
  * Schmidt A59 (2021a) Schmidt R. E., 2021a, CBET, 4906
  * Schmidt A60 (2021b) Schmidt R. E., 2021b, JAAVSO, 49, 1
  * Schreiber et al. (2016) Schreiber M. R., Zorotovic M., Wijnen T. P. G., 2016, MNRAS, 455, L16
  * Shafter A61 (2006) Shafter A. W., Misselt K. A., 2006, ApJ, 644, 1104
  * Shafter A63 (1993) Shafter A. W., Misselt K. A., Veal J. M., 1993, PASP, 105, 853
  * Shafter A62 (1995) Shafter A. W., Misselt K. A., Szkody P., Politano M., 1995, ApJ, 448, L33
  * Shara A64 (2007) Shara M. M. et al. 2007, Nature, 446, 159
  * Shara A65 (2012) Shara M. M. et al., 2012, ApJ, 756, 107
  * Shara A66 (2017a) Shara M. M. et al., 2017a, Nature, 548, 558
  * Shara A67 (2017b) Shara M. M., Drissen L., Martin T., Alarie A., Stephenson F. R., 2017b, MNRAS, 465, 739
  * Shara et al. (2018) Shara M. M., Prialnik D., Hillman Y., Kovetz A., 2018, ApJ, 860, 110
  * Soszynski A68 (2016) Soszynski I. et al., 2016, Acta Astronomica, 66, 405
  * Strope, Schaefer, and Henden (2010) Strope R. J., Schaefer B. E., Henden A. A., 2010, AJ, 140, 34
  * Surina et al. (2011) Surina F., Bode M. F., Darnley M. J., 2011, arXiv:1111.5524
  * Swierczynski A69 (2010) Swierczynski E., et al., 2010, in Prsa A., Zejda M., eds, Astronomical Society of the Pacific Conference Series Vol. 435, Binaries - Key to Comprehension of the Universe. p. 297 (arXiv:0909.5149)
  * Szkody (1994) Szkody P., 1994, AJ, 108, 639
  * Tappert et al. (2013) Tappert C., Schmidtobreick L., Vogt N., Ederoclite A., 2013, MNRAS, 436, 2412
  * Thomas A70 (2021) Thomas N., Ziegler K., Liu P., 2021, JAAVSO, 49, 1
  * Thorstensen A72 (1995) Thorstensen J. R., Ringwald F. A., 1995, IBVS, 4249
  * Thorstensen A74 (2000) Thorstensen J. R., Taylor C. J., 2000, MNRAS, 312, 629
  * Thorstensen A71 (2010) Thorstensen J. R., Peters C. S., Skinner J. N., 2010, PASP, 122, 1285
  * Thorstensen A73 (2017) Thorstensen J. R., Ringwald F. A., Taylor C. J., Sheets H. A., Peters C. S., Skinner J. N., Alper E. H., Weil K. E., 2017, RNAAS, 1, 29
  * van Belle et al. (1999) van Belle G. T. et al., 1999, AJ, 117, 521
  * Vogt et al. (2018) Vogt N., Tappert C., Puebla E. C., Fuentes-Morales I., Ederoclite A., Schmidtobreick L., 2018, MNRAS, 478, 5427
  * Wagner A75 (2012) Wagner R. M. et al., 2012, ATel, 4157
  * Walker (1954) Walker M. F., 1954, PASP, 66, 230
  * Walter A76 (2015) Walter F. M., 2015, in Dufour P., Bergeron P., Fontaine G., eds, Astronomical Society of the Pacific Conference Series Vol. 493, 19th European Workshop on White Dwarfs. p. 507
  * Walter et al. (2012) Walter F. M., Battisti A., Towers S. E., Bond H. E., Stringfellow G. S., 2012, PASP, 124, 1057
  * Warner A77 (2009) Warner B., Woudt P. A., 2009, MNRAS, 397, 979
  * Weight A78 (1994) Weight A., Evans A., Naylor T., Wood J., Bode M., 1994, MNRAS, 266, 761
  * Williams et al. (1991) Williams R. E., Hamuy M., Phillips M. M., Haethcote S. R., Wells L., Navarrete M., 1991, ApJ, 376, 721
  * Williams et al. (2016) Williams S. C., Darnley M. J., Bode M. F., Shafter A. W., 2016, ApJ, 817, 143
  * Woudt A79 (2001) Woudt P. A., Warner B., 2001, MNRAS, 328, 159
  * Woudt and Warner (2002) Woudt P. A., Warner B., 2002, MNRAS, 335, 44
  * Woudt A80 (2003) Woudt P. A., Warner B., 2003, MNRAS, 340, 1011
  * Woudt et al. (2005) Woudt P. A., Warner B., Spark M., 2005, MNRAS, 364, 107
  * Woudt et al. (2009) Woudt P. A., Warner B., Osborne J., Page K., 2009, MNRAS, 395, 2177
  * Zemko A81 (2018) Zemko P. et al., 2018, MNRAS, 480, 4489