#  On the Schur multiplier of nilpotent Lie superalgebra

Afsaneh Shamsaki  School of Mathematics and Computer Science
Damghan University, Damghan, Iran<EMAIL_ADDRESS>and  Peyman
Niroomand<EMAIL_ADDRESS><EMAIL_ADDRESS>School of Mathematics
and Computer Science
Damghan University, Damghan, Iran

###### Abstract.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra where $m+n\geq 4$
and $n\geq 1$. This paper classifies such nilpotent Lie superalgebras $L$ with
a derived subsuperalgebra of dimension $m+n-2$ such that
$\gamma(L)=m+2n-2-\dim\mathcal{M}(L)$, where $\gamma(L)\in\\{0,1,2\\}$ and
$\mathcal{M}(L)$ denotes the Schur multiplier of $L$. Furthermore, we show
that all these superalgebras are capable.

###### Key words and phrases:

Schur multiplier; nilpotent Lie superalgebra

Mathematics Subject Classification 2010. 17B30

##  1\. Introduction and Preliminaries

Let $L$ be a finite-dimensional Lie superalgebra, and let $L\cong F/R$, where
$F$ is a free Lie superalgebra. The Schur multiplier $\mathcal{M}(L)$ of $L$
is isomorphic to $R\cap F^{2}/[R,F]$ as shown in [3]. Classifying Lie
superalgebras of arbitrary dimensions without any constraints is a challenging
task. In fact, the classification of all finite-dimensional Lie superalgebras
remains an open problem. A common approach is to classify nilpotent Lie
superalgebras based on a fixed invariant. Some researchers focus on the
structure of $L$ by examining the dimension of the Schur multiplier. This
invariant aids in classifying certain nilpotent Lie superalgebras. For
example, the classification of nilpotent Lie superalgebras based on the
dimension of the Schur multiplier is discussed in [4]. According to [3], the
dimension of the Schur multiplier of an $(m\mid n)$-dimensional nilpotent Lie
superalgebra $L$ with the derived supersubalgebra of maximum dimension $m+n-2$
is bounded by $m+2n-2$. Thus, there exists a non-negative integer $\gamma(L)$
such that $\gamma(L)=m+2n-2-\dim\mathcal{M}(L)$.

Throughout this paper, we assume that $L$ is a Lie superalgebra over a field
of characteristic not equal to 2 or 3. We classify all finite-dimensional
nilpotent Lie superalgebras $L$ when $\gamma(L)\in\\{0,1,2\\}$. Additionally,
we demonstrate that all these Lie superalgebras are capable. Recall from [7]
that a Lie superalgebra $L$ is said to be capable if there exists a Lie
superalgebra $M$ such that $L\cong M/Z(M)$. Also, $Z^{*}(L)$ is the epicenter
of a Lie superalgebra $L$, defined as the smallest ideal in $L$ such that
$L/Z^{*}(L)$ is capable. Clearly, the ideal $Z^{*}(L)=0$ if and only if $L$ is
capable. For more details on capable Lie superalgebras, see [7].

The following theorem is a useful tool for identifying capable Lie
superalgebras.

######  Theorem 1.1.

[6, Theorem 4.9] Let $N$ be a central ideal of the Lie superalgebra $L$. Then
the following conditions are equivalent:

  * (i)

$\mathcal{M}(L)\cong\mathcal{M}(L/N)/N\cap L^{2},$

  * (ii)

$N\subseteq Z^{*}(L),$

  * (iii)

the natural map $\mathcal{M}(L)\longrightarrow\mathcal{M}(L/N)$ is a
monomorphism.

The lower central series is defined by $C^{0}(L)=L$ and inductively by
$C^{i+1}(L)=[L,C^{i}(L)]$. A Lie superalgebra $L$ is nilpotent if there exists
an integer $n$ such that $C^{n}(L)=0$. Let $L=L_{\overline{0}}\oplus
L_{\overline{1}}$. Similarly, we may define the following sequences for $L$ as
defined in [2]:

$\displaystyle C^{0}(L_{\overline{0}})=L_{\overline{0}}\quad\text{and}\quad
C^{k}(L_{\overline{0}})=[L_{\overline{0}},C^{k-1}(L_{\overline{0}})]\quad\text{for
all}\quad k\geq 1,$ $\displaystyle
C^{0}(L_{\overline{1}})=L_{\overline{1}}\quad\text{and}\quad
C^{k}(L_{\overline{1}})=[L_{\overline{0}},C^{k-1}(L_{\overline{1}})]=ad_{L_{\overline{0}}}(C^{k-1}(L_{\overline{1}}))\quad\text{for
all}\quad k\geq 1.$

######  Theorem 1.2.

[2, Theorem 2.1] The Lie superalgebra $L=L_{\overline{0}}\oplus
L_{\overline{1}}$ is nilpotent if and only if there exist
$(p,q)\in\mathbb{N}\times\mathbb{N}$ such that $C^{p-1}(L_{\overline{0}})\neq
0$, $C^{q-1}(L_{\overline{1}})\neq 0$, and
$C^{p}(L_{\overline{0}})=C^{q}(L_{\overline{1}})=0$.

For convenience, we list some results regarding the Schur multiplier of a Lie
superalgebra.

######  Theorem 1.3.

[3, Theorem 3.4] Let $L$ be an $(m|n)$-dimensional Lie superalgebra. Then

$\dim\mathcal{M}(L)=\frac{1}{2}\big{(}(m+n)^{2}+(n-m)\big{)}$

if and only if $L$ is abelian.

######  Theorem 1.4.

[3, Theorem 3.7] Let $L=L_{\overline{0}}\oplus L_{\overline{1}}$ be a Lie
superalgebra with $\dim L=(m|n)$ and let $K\subseteq Z(L)$ be a graded ideal.
If $H\cong L/K$ is the quotient Lie superalgebra, then

$\displaystyle(i)\quad\dim\mathcal{M}(L)+\dim(L^{2}\cap
K)\leq\dim\mathcal{M}(H)+\dim\mathcal{M}(K)+\dim(H/H^{2}\otimes K/K^{2}),$
$\displaystyle(ii)\quad\dim\mathcal{M}(L)+\dim(L^{2}\cap
K)\leq\frac{1}{2}\big{(}(m+n)^{2}+(n-m)\big{)}.$

######  Theorem 1.5.

[3, Theorem 5.4] Let $L$ be a non-abelian nilpotent Lie superalgebra with
$\dim L=(m|n)$, where $m+n\geq 3$ and $\dim L^{2}=m+n-2$.

  * (i)

If $m+n=3$, then $\dim\mathcal{M}(L)$ is equal to 1 or 2.

  * (ii)

If $m+n\geq 4$ and $n\geq 1$, then $\dim\mathcal{M}(L)\leq m+2n-2$.

##  2\. Main Result

In this section, we classify the structure of all $(m|n)$-dimensional
nilpotent Lie superalgebras with a derived subsuperalgebra of dimension
$m+n-2$, for $m+n\geq 4$, $n\geq 1$, and $\gamma(L)\in{0,1,2}$. Furthermore,
we demonstrate that all such Lie superalgebras are capable.

The following proposition provides the Schur multipliers of
$(m|n)$-dimensional nilpotent Lie superalgebras with a derived subsuperalgebra
of dimension $m+n-2$ for $4\leq m+n\leq 5$ and $n\geq 1$.

######  Proposition 2.1.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra with derived
subalgebra of dimension $m+n-2$ such that $4\leq m+n\leq 5$ and $n\geq 1.$
Then we have the following table.

Table 1.  Name | Basis | Nonzero multiplication | $\dim\mathcal{M}(L)$
---|---|---|---
$(2|2)_{1}$ | $e_{1},e_{2},f_{1},f_{2}$ | $[f_{1},f_{1}]=e_{1},[f_{2},f_{2}]=e_{2}$ | $1$
$(2|2)_{4}$ | $e_{1},e_{2},f_{1},f_{2}$ | $[f_{1},f_{2}]=e_{1},[f_{2},f_{2}]=e_{2}$ | $2$
$(2|2)_{6}$ | $e_{1},e_{2},f_{1},f_{2}$ | $[e_{2},f_{2}]=f_{1},[f_{2},f_{2}]=e_{1}$ | $2$
$(1|3)_{1}$ | $e_{1},f_{1},f_{2},f_{3}$ | $[e_{1},f_{2}]=f_{1},[e_{1},f_{3}]=f_{2}$ | $3$
$(1|4)_{7}$ | $e_{1},f_{1},...,f_{4}$ | $[e_{1},f_{2}]=f_{1},[e_{1},f_{3}]=f_{2},[e_{1},f_{4}]=f_{3}$ | $3$
$(3|2)_{5}$ | $e_{1},e_{2},e_{3},f_{1},f_{2}$ | $[f_{1},f_{1}]=e_{2},[f_{1},f_{2}]=e_{1},[f_{2},f_{2}]=e_{3}$ | $2$
$(3|2)_{13}$ | $e_{1},e_{2},e_{3},f_{1},f_{2}$ | $[e_{1},e_{2}]=e_{3},[e_{1},f_{2}]=f_{1},[f_{1},f_{2}]=e_{3}$ | $3$
|  | $[f_{2},f_{2}]=2e_{2}$ |
$(2|3)_{18}$ | $e_{1},e_{2},f_{1},f_{2},f_{3}$ | $[e_{1},f_{3}]=f_{1},[e_{2},f_{2}]=f_{1},[f_{2},f_{2}]=2e_{1}$ | $2$
|  | $[f_{2},f_{3}]=-e_{2}$ |
$(2|3)_{19}$ | $e_{1},e_{2},f_{1},f_{2},f_{3}$ | $[e_{1},f_{3}]=f_{1},[e_{2},f_{2}]=f_{1},[f_{2},f_{3}]=-e_{1}$ | $3$
|  | $[f_{3},f_{3}]=2e_{2}$ |
$(2|3)_{22}$ | $e_{1},e_{2},f_{1},f_{2},f_{3}$ | $[e_{1},f_{2}]=f_{1},[e_{1},f_{3}]=f_{2},[f_{3},f_{3}]=e_{2}$ | $3$
$(2|3)_{23}$ | $e_{1},e_{2},f_{1},f_{2},f_{3}$ | $[e_{1},f_{2}]=f_{1},[e_{1},f_{3}]=f_{2},[f_{1},f_{3}]=-e_{2}$ | $2$
|  | $[f_{2},f_{2}]=e_{2}$ |

###### Proof.

Let $L\cong(2|3)_{22}=\langle
e_{1},e_{2},f_{1},f_{2},f_{3}\mid[e_{1},f_{2}]=f_{1},[e_{1},f_{3}]=f_{2},[f_{3},f_{3}]=e_{2}\rangle$.
By using the method of in [3], we compute $\dim\mathcal{M}((2|3)_{22})$. Start
with

$\displaystyle[e_{1},e_{2}]=s_{1},$ $\displaystyle[e_{1},f_{1}]=s_{2},$
$\displaystyle[e_{1},f_{2}]=f_{1}+s_{3},$
$\displaystyle[e_{1},f_{3}]=f_{2}+s_{4},$ $\displaystyle[e_{2},f_{1}]=s_{5},$
$\displaystyle[e_{2},f_{1}]=s_{6},$ $\displaystyle[e_{2},f_{3}]=s_{7},$
$\displaystyle[f_{1},f_{2}]=s_{8},$ $\displaystyle[f_{1},f_{3}]=s_{9},$
$\displaystyle[f_{1},f_{1}]=s_{10},$ $\displaystyle[f_{2},f_{3}]=s_{11},$
$\displaystyle[f_{2},f_{2}]=s_{12},$
$\displaystyle[f_{3},f_{3}]=e_{2}+s_{13},$

where the set $\\{s_{1},...,s_{13}\\}$ generates $\mathcal{M}(L)$. The graded
Jacobi identity on all possible triples gives

$\displaystyle
s_{1}=[e_{1},[f_{3},f_{3}]]=[f_{3},[e_{1},f_{3}]]-[f_{3},[f_{3},e_{1}]]=2s_{11},$
$\displaystyle
s_{5}=[e_{2},[e_{1},f_{2}]]=-[f_{2},[e_{2},e_{1}]]-[e_{1},[f_{2},e_{2}]]=0,$
$\displaystyle
s_{6}=[f_{3},[e_{2},e_{1}]]=[e_{1},[f_{3},e_{2}]]+[x_{2},[x_{5},x_{1}]]=0,$
$\displaystyle
s_{7}=[f_{3},[f_{3},f_{3}]]=-[f_{3},[f_{3},f_{3}]]-[f_{3},[f_{3},f_{3}]]=-2s_{7},$
$\displaystyle
s_{8}=[f_{1},[e_{1},f_{3}]]=[f_{3},[f_{1},e_{1}]]+[e_{1},[f_{3},f_{1}]]=0,$
$\displaystyle
s_{9}=[f_{1},[e_{1},f_{3}]]=[f_{3},[f_{1},e_{1}]]+[e_{1},[f_{3},f_{1}]]=0,$
$\displaystyle
s_{10}=[f_{1},[e_{1},f_{2}]]=[f_{2},[f_{1},e_{1}]]+[e_{1},[f_{2},f_{1}]]=0.$

Let $f^{\prime}{1}=f{1}+s_{3}$, $f^{\prime}{2}=f{2}+s_{4}$, and
$e^{\prime}{2}=e{2}+s_{13}$. A change of variables allows
$s_{3}=s_{4}=s_{13}=0$. Hence $\mathcal{M}((2|3){22})=\langle
s{1},s_{2},s_{9}\rangle$, and thus $\dim\mathcal{M}((2|3)_{22})=3$.

By a similar method, the Schur multiplier of other Lie superalgebras in Table
LABEL:tab1 is obtained. ∎

In the remainder of the paper, we consider all $(m|n)$-dimensional nilpotent
Lie superalgebras $L$ with the derived subsuperalgebra of dimension $m+n-2$
such that $m+n\geq 4$ and $n\geq 1.$

Here’s the revised draft:

“‘latex

######  Lemma 2.2.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra with the Schur
multiplier of dimension $m+2n-2-i$ for all $i\geq 0$, and let $K\subseteq
L^{2}\cap Z(L)$ be an $(1,0)$-dimensional ideal. Then $\dim(L/K)^{2}=m+n-3$,
and

  * (i)

if $m+n=4$, then $\dim\mathcal{M}(L/K)$ is equal to $1$ or $2$,

  * (ii)

if $m+n\geq 5$, then $m+2n-2-(i+1)\leq\dim\mathcal{M}(L/K)\leq m+2n-3$ for all
$i\geq 0$.

###### Proof.

Assume that $\dim L^{2}=m+n-2$ and $\dim L^{2}\cap K=1$. Thus,
$\dim(L/K)^{2}=m+n-3$.
(i) Since $\dim K=(1,0)$, we have $\dim L/K=(m-1,n)$ such that $m+n=3$, and by
part (i) of Theorem 1.5, $\dim\mathcal{M}(L/K)=1$ or $2$.
(ii) We know that $\dim\mathcal{M}(K)=0$ by Theorem 1.3. It implies that
$m+2n-(i+1)\leq\dim\mathcal{M}(L/K)$ for all $i\geq 0$, by Theorem 1.4. On the
other hand, $\dim\mathcal{M}(L/K)\leq m+2n-3$ by part (ii) of Theorem 1.5.
Hence the result follows. ∎

######  Lemma 2.3.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra with the Schur
multiplier of dimension $m+2n-2-i$ for all $i\geq 0$, and let $K\subseteq
L^{2}\cap Z(L)$ be an $(0,1)$-dimensional ideal. Then $\dim(L/K)^{2}=m+n-3$,
and

  * (i)

if $m+n=4$, then $\dim\mathcal{M}(L/K)$ is equal to $1$ or $2$,

  * (ii)

if $m+n\geq 5$ and $n=1$, then $\dim\mathcal{M}(L/K)\leq m-2$,

  * (iii)

if $m+n\geq 5$ and $n\geq 2$, then $m+2n-2-(i+2)\leq\dim\mathcal{M}(L/K)\leq
m+2n-4$ for all $i\geq 0$.

###### Proof.

The parts (i) and (iii) are obtained by a similar method to Lemma 2.2. It is
sufficient to prove part (ii).
(ii) Since $m+n\geq 5$ and $n=1$, we have $\dim(L/K)=(m,0)$. So, $L/K$ is a
Lie algebra, and $\dim\mathcal{M}(L/K)\leq m-2$ by [5, Theorem 3.1]. ∎

######  Theorem 2.4.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra. Then
$\dim\mathcal{M}(L)<m+2n-3$.

###### Proof.

We proceed by induction on $m+n$. Let $m+n=4$. By looking at Table LABEL:tab1,
there is no $(m|n)$-Lie superalgebra such that $m+n=4$ and $n\geq 1$ with
$\dim\mathcal{M}(L)<m+2n-2$. Let $m+n\geq 5$ and $K$ be a central ideal of $L$
such that $\dim K=1$ and $K\subseteq Z(L)\cap L^{2}$. By Theorem 1.4, we have

$\displaystyle
1+\dim(\mathcal{M}(L))\leq\dim(\mathcal{M}(L/K))+\dim(\mathcal{M}(K))+\dim(L/L^{2}\otimes
K)$

If $K\subseteq L_{\overline{0}}$, then $\dim(\mathcal{M}(K))=0$ by Theorem
1.3. Also, we have $\dim(\mathcal{M}(L/K))<m+2n-3$ by the induction
hypothesis. So, $\dim\mathcal{M}(L)<(m+2n-3)-1+\dim L/L^{2}$ by Theorem 1.4.
If $K\subseteq L_{\overline{1}}$, then $\dim\mathcal{M}(K)=1$ by Theorem 1.3.
Hence

$\displaystyle\dim\mathcal{M}(L)<(m+2n-4)+\dim L/L^{2}$

Thus $\dim\mathcal{M}(L)<m+2n-2$. Suppose that $\dim\mathcal{M}(L)=m+2n-3$. By
a similar method, we can see $\dim\mathcal{M}(L)<m+2n-3$. The result follows.
∎

Here’s the revised version:

######  Proposition 2.5.

No $(m|n)$-dimensional nilpotent Lie superalgebras exist with a derived
subsuperalgebra of dimension $4$, and $\dim\mathcal{M}(L)=m+2n-4$ when $m+n=6$
and $n\geq 1$.

###### Proof.

Given $m+n=6$ and $n\geq 1$, the superalgebra $(m|n)$ can be one of the
following: $(5|1)$, $(4|2)$, $(3|3)$, $(2|4)$, $(1|5)$, or $(0|6)$. Assuming
$K\subseteq L^{2}\cap Z(L)$ and $(m|n)=(4|2)$, if $K\subseteq
L_{\overline{0}}$, then $\dim L/K=(3|2)$. By Lemma 2.2, $\dim(L/K)^{2}=3$, and
$\dim\mathcal{M}(L/K)\geq 3$. Thus, $L/K$ is isomorphic to a nilpotent Lie
superalgebra of superdimension $(3|2)_{13}$ as per Table LABEL:tab1. Let
$K=\langle e_{4}\rangle$ and

$(3|2)_{13}=\langle
e_{1},e_{2},e_{3},f_{1},f_{2}|[e_{1},e_{2}]=e_{3},[e_{1},f_{2}]=f_{1},[f_{1},f_{2}]=e_{3},[f_{2},f_{2}]=2e_{2}\rangle.$

We obtain $L=\langle e_{1},e_{2},e_{3},e_{4},f_{1},f_{2}\rangle$ with
brackets:

$\displaystyle[e_{1},e_{2}]=e_{3}+\alpha_{1}e_{4},~{}\quad[e_{1},e_{3}]=\alpha_{2}e_{4},\quad[e_{1},f_{1}]=0,$
$\displaystyle[e_{1},f_{2}]=f_{1},~{}\qquad\qquad[e_{2},e_{3}]=\alpha_{3}e_{4},\quad[e_{2},f_{1}]=0,$
$\displaystyle[e_{2},f_{2}]=0,\qquad\quad\qquad[e_{3},f_{1}]=0,~{}\qquad[e_{3},f_{2}]=0,$
$\displaystyle[f_{1},f_{2}]=e_{3}+\alpha_{4}e_{4},~{}\quad[f_{1},f_{1}]=\alpha_{5}e_{4},\quad[f_{2},f_{2}]=2e_{2}+\alpha_{6}e_{4}.$

Let $e^{\prime}_{3}=e_{3}+\alpha_{1}e_{4}$ and
$e^{\prime}_{2}=e_{2}+1/2\alpha_{6}e_{4}$. A change of variables sets
$\alpha_{1}=\alpha_{6}=0$. By graded Jacobi identities, we have

$\displaystyle\alpha_{2}e_{4}=[e_{1},e_{3}]=[e_{1},[f_{1},f_{2}]]=[f_{2},[e_{1},f_{1}]]-[f_{1},[f_{2},e_{1}]]=\alpha_{5}e_{4},$
$\displaystyle\alpha_{3}e_{4}=[e_{2},e_{3}]=[e_{2},[f_{1},f_{2}]]=[f_{2},[e_{2},f_{1}]]-[f_{1},[f_{2},e_{2}]]=0.$

Thus,

$L=\langle
e_{1},e_{2},e_{3},e_{4},f_{1},f_{2}|[e_{1},e_{2}]=e_{3},[e_{1},e_{3}]=[f_{1},f_{1}]=\alpha_{2}e_{4},[e_{1},f_{2}]=f_{1},$
$[f_{1},f_{2}]=e_{3}+\alpha_{4}e_{4},[f_{2},f_{2}]=2e_{2}\rangle.$

Since $\dim L^{2}=4$, $\alpha_{2}$ or $\alpha_{4}$ must not be zero. By a
similar method as in Proposition 2.1, $\dim\mathcal{M}(L)\leq 3$. Thus,
there’s no such Lie superalgebra of dimension $(4|2)$ with
$\dim\mathcal{M}(L)=4$. If $K\subseteq L_{\overline{1}}$, then $\dim
L/K=(4|1)$. By Table LABEL:tab1, there’s no such Lie algebra.

For $(m|n)=(5|1)$, $(3|3)$, $(2|4)$, $(1|5)$, or $(0|6)$, there’s no such
nilpotent Lie superalgebra $L$ with derived subsuperalgebra of dimension $4$
such that $m+n=6$, $n\geq 1$, and $\dim\mathcal{M}(L)=m+2n-4$ by a similar
method. ∎

######  Theorem 2.6.

For an $(m|n)$-dimensional nilpotent Lie superalgebra with $m+n\geq 6$ and
$n\geq 1$, $\dim\mathcal{M}(L)<m+2n-4$.

###### Proof.

The proof follows similar methods as in the proof of Theorem 2.4, Proposition
2.5, Theorems LABEL:tab1 and 1.4. ∎

We now derive the structure of all nilpotent Lie superalgebras of dimension
$(m|n)$ when $\gamma(L)\in\\{0,1,2\\}$.

######  Main Theorem.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra.

  * (i)

There is no such Lie superalgebra with $\gamma(L)=0.$

  * (ii)

There is no such Lie superalgebra with $\gamma(L)=1.$

  * (iii)

$\gamma(L)=2$ if and only if $L$ is isomorphic to one of the Lie superalgebras
$(2|2)_{4},$ $(2|2)_{6},$ $(1|3)_{1},$ $(3|2)_{13}$ or $(2|3)_{18}.$

###### Proof.

Suppose $\gamma(L)=0$ or $1.$ Then $\dim\mathcal{M}(L)=m+2n-2$ or $m+2n-3,$
respectively. Therefore, by Theorem 2.4, there exists no $(m|n)$-dimensional
nilpotent Lie superalgebra with $4\leq m+n\leq 5$ and $n\geq 1$.

For part (iii), using Theorem 2.6, there is no such nilpotent Lie superalgebra
$L$ with $\dim\mathcal{M}(L)=m+2n-4$ and $\dim L^{2}=m+n-2$ for $m+n\geq 6.$
Hence, for every nilpotent Lie superalgebra with $\dim\mathcal{M}(L)=m+2n-4$
and $\dim L^{2}=m+n-2$, we have $m+n\leq 5.$ By examining Table LABEL:tab1,
$L$ should be isomorphic to one of the Lie superalgebras $(2|2)_{4},$
$(2|2)_{6},$ $(1|3)_{1},$ $(3|2)_{13}$ or $(2|3)_{18}.$ The converse is
obtained by referring to Table LABEL:tab1. ∎

######  Corollary 2.7.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra with $n+m\geq 6$.
Then

$\dim\mathcal{M}(L)\leq m+2n-5.$

######  Corollary 2.8.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra. If
$\gamma(L)=2,$ then $L$ is capable.

###### Proof.

Let $L$ be an $(m|n)$-dimensional nilpotent Lie superalgebra with the derived
supersubalgebra of dimension $m+n-2$ and $\dim\mathcal{M}(L)=m+2n-4.$ Then
$\dim\mathcal{M}(L/K)\leq\dim\mathcal{M}(L)$ by Lemmas 2.2 and 2.3. Hence,
$(2|2)_{4},$ $(2|2)_{6},$ $(1|3)_{1},$ $(3|2)_{13}$ and $(2|3)_{18}$ are
capable by Theorem 1.1. The result follows. ∎

## References

  * [1] Alvarez, M. A. and Hernández, I.; Varieties of nilpotent Lie superalgebras of dimension $\leq 5$. Forum Math. 32 (2020), no. 3, 641–661.
  * [2] Gilg, M.; Low-dimensional filiform Lie superalgebras. Revista Matemática Complutense, 14(2), (2001), 463–478.
  * [3] Nayak, S.; Multipliers of nilpotent Lie superalgebras. Comm. Algebra 47 (2019), no. 2, 689–705.
  * [4] Nayak, S.; Classification of finite dimensional nilpotent Lie superalgebras by their multipliers. J. Lie Theory 31.2 (2021), 439–458.
  * [5] Niroomand, P. and Russo, F.; A restriction on the Schur multiplier of nilpotent Lie algebras. The Electronic Journal of Linear Algebra 22 (2011), 1–9.
  * [6] Padhan, R. N. and Nayak, S.; On capability and the Schur multiplier of some nilpotent Lie superalgebras. Linear and Multilinear Algebra (2020), 1–12.
  * [7] Padhan, R. N., Nayak, S. and Pati, K. C.; Detecting capable Lie superalgebras. Communications in Algebra 49.10, (2021), 4274–4290.


# Improving MATLAB’s _isprime_ performance without arbitrary-precision
arithmetic

Travis Near
Department of Computer Science
Missouri University of Science and Technology
Rolla, MO 65409
<EMAIL_ADDRESS>

(August 1, 2021)

###### Abstract

MATLAB is a numerical computing platform used by scientists, engineers,
mathematicians, and students which contains many mathematical functions,
including _isprime_. MATLAB’s _isprime_ function determines which elements of
an input array are prime. This research details modular arithmetic techniques,
the Miller—Rabin primality test, vectorized operations, and division-
minimizing strategies which harness the power of MATLAB’s capabilities to
improve _isprime_ ’s performance. The results are typically 5 to 10 times
faster for small integers and many hundreds of times faster for large integers
and long arrays.

_Keywords_ MATLAB $\cdot$ _isprime_ $\cdot$ _isprime_fast_ $\cdot$ performance

##  1 Introduction

MATLAB’s _isprime_ function accepts an input array, _N_ , then returns which
elements in the array are prime:

        >> isprime([1, 2, 3, 4, 5])
        ans =

          1x10 logical array

          0   1   1   0   1

_isprime_ begins with the sieve of Eratosthenes followed by trial division of
all primes up to the square root of the maximum value, _M_ , in the array [5].
The sieve begins with an incrementing sequence of $[2,3,4,5,...,\emph{M}]$.
Then it _crosses-off_ multiples of primes. Two, for example, is prime so all
multiples of two are crossed-off given that they share a common factor. This
leaves $[2,3,5,7,9,...,\emph{M}]$. The next prime is three: it is marked prime
then its multiples are crossed-off leaving $[2,3,5,7,11,...,\emph{M}]$. The
process repeats until all prime numbers have been visited.

Next, _isprime_ utilizes trial division of all primes returned by the sieve.
_isprime_ divides each element of the input array by all primes which are
_less_ than that element. The input element is marked prime when none of the
smaller primes divides into it evenly. This algorithm is simple, effective,
and accurate, but is not optimized for any category of inputs.

This research introduces _isprime_fast_ , a replacement for _isprime_.
_isprime_fast_ ’s inputs and outputs are identical; the only difference lies
in the techniques used to determine primality. _isprime_fast_ specializes for
four classifications of inputs: (1) large scalar integers, (2) medium scalar
integers, (3) small scalar integers, and (4) arrays.

##  2 Improving the performance for large scalar integers (249 to 264)

MATLAB’s _isprime_ takes about 30 seconds to check primality for 64-bit
integers. _isprime_fast_ reduces this to about 0.1 seconds by implementing a
modified Miller–Rabin primality test. MATLAB does not support arbitrary
precision arithmetic, so the modifications made to Miller–Rabin use modular
arithmetic properties to clamp the range of intermediate values to ensure a
correct result.

Miller–Rabin is a _probabilistic_ prime number algorithm, yet _isprime_ is
_deterministic_. To replace _isprime_ , _isprime_fast_ must use a
deterministic version of Miller–Rabin. The reason that Miller–Rabin is
probabilistic is due to its selection of _bases_ to check congruence relations
[10]. To produce deterministic results, Miller–Rabin needs to select bases
which produce the correct result for _isprime_ ’s entire input domain of 0 to
264.

###  2.1 The Miller–Rabin primality test

Miller–Rabin is built on two properties which hold for all prime numbers.
First, Fermat’s little theorem [10] which states that for an integer base $a$
and prime number $n$,

$a^{n-1}\equiv\pm 1\pmod{n}$ (1)

Second, the square roots of Equation 1 also hold. By rewriting $n-1$ as powers
of two multiplied by a constant $d$, this becomes,

$a^{2^{s}d}\equiv\pm 1$ (2)

By successively removing powers of two, the terms to check for congruence are,

$a^{2^{s}d},a^{2^{s-1}d},a^{2^{s-2}d},...,a^{d}$ (3)

If $n$ is prime, all elements in this sequence are congruent to $\pm 1$ (mod
$n$) [4].

These two properties hold for all prime numbers. However, there exist
composite numbers that pass these two tests which are known as _strong
pseudoprimes_. The likelihood that a number is truly prime increases as more
bases are checked against the two properties above. By carefully selecting
bases, Miller–Rabin can be made deterministic up to any arbitrary value. To
replace _isprime_ , Miller–Rabin must use bases which return the correct
answer for inputs up to 264. The following bases provide accurate results for
all 64-bit integers [13]:

$2,325,9375,28178,450775,9780504,1795265022$ (4)

_isprime_fast_ implements the Miller–Rabin algorithm as detailed by its
original authors. However, it makes several modifications to ensure that no
intermediate result at any point exceeds MATLAB’s maximum integer value of
264. This adapted Miller–Rabin algorithm is built on three core functions:
_ModAdd_ , _ModMultiply_ , and _ModExp_.

###  2.2 ModAdd

_ModAdd_ calculates $c=a+b\pmod{m}$. This function adds two unsigned 64-bit
integers followed by (_mod m_) without overflow. The temporary result of $a+b$
(before _mod m_) can overflow (more specifically, _saturate_) and return the
highest possible value given the precision available. For example,
$100+200=300$ but for unsigned 8-bit integers (_uint8_) the result is the
largest value allowed by this data type: $100+200=255$ for _uint8_ in MATLAB.
Saturation is undesirable for _isprime_fast_ because it surrenders accuracy.
To avoid overflowing, _ModAdd_ uses modular arithmetic properties, invariants,
and careful ordering of the sequence of its calculations:

  * •

$c=a+b\pmod{m}$ _[intended calculation]_

  * •

$assert(a+b>m)$

  * •

$assert(m\geq b)$

  * •

$assert(m\geq a)$

  * •

$assert(m\leq 2^{64})$

  * •

Given that $a+b>m$, the inequality $a>m-b$ must hold

  * •

Given that $m\geq b$, the result $m-b$ must be positive or zero

  * •

$c=a+b-m$ _[mod can be simplified to subtraction given that a and b are less
than m]_

  * •

$c=a-(m-b)$ _[Property 1]_

###  2.3 ModMultiply

_ModMultiply_ calculates $c=a*b\pmod{m}$. This function multiplies two
unsigned 64-bit integers followed by _mod m_ without overflow. Like _ModAdd_ ,
the temporary result of $a*b$ before _mod m_ will overflow when $a*b$ exceeds
the range of 64-bit integers. An optimized implementation of _ModMultiply_
considers four scenarios:

  1. 1.

Fast path: $a*b$ does not overflow

    1. (a)

$c=a*b\pmod{m}$ _[intended calculation]_

  2. 2.

$a=b$

    1. (a)

$c=a*b\pmod{m}$

    2. (b)

$assert(m\geq a)$

    3. (c)

$c=a*a\pmod{m}$ _[substitute a=b]_

    4. (d)

$c=(a\pmod{m}*a\pmod{m})\pmod{m}$ _[modular distribution]_

    5. (e)

$c=(m-a)*(m-a)\pmod{m}$ _[mod is based on difference between m and a]_

    6. (f)

$c=(a*a)\pmod{m}=(m-a)*(m-a)\pmod{m}$

    7. (g)

$aSmall=min(m,m-a)$ _[from line f, m and m – a give the same result]_

    8. (h)

if $aSmall\leq 2^{32}$; $c$ can be calculated without overflow using _uint64_
, $c=aSmall^{2}\pmod{m}$

    9. (i)

otherwise, proceed to case 4 (general overflow)

  3. 3.

$a=2$

    1. (a)

$c=a*b\pmod{m}$ _[intended calculation]_

    2. (b)

$assert(m\geq b)$

    3. (c)

$assert(2*b>m)$ _[must be true given that $a*b>2^{64}>m$]_

    4. (d)

$c=2*b\pmod{m}$ _[substitute a=2]_

    5. (e)

$c=b+b\pmod{m}$ _[replace multiplication with addition]_

    6. (f)

$c=b-(m-b)$ _[by Property 1]_

  4. 4.

$a*b$ overflows 64-bit integers

    1. (a)

A common modular-multiply implementation [2] (shown below as
_ModMultiplyClassic_) computes $a*b\pmod{m}$ by reducing to
$mod(2*a,m)*(b/2)$.

    2. (b)

This increases _ModMultiply_ ’s range, but still encounters an overflow when
$2^{63}<a\leq 2^{64}$. This is due to doubling $a$ before applying _mod_.

    3. (c)

To replace _isprime_ , _ModMultiply_ ’s range must be doubled to reach
$2^{64}(2*2^{63})$.

⬇

1function c = ModMultiplyClassic(a, b, m)

2 % c = a * b (mod m)

3 c = 0;

4 a = mod(a, m);

5 while b > 0

6 if bitand(b, 1) % odd numbers

7 c = mod(c + a, m);

8 end

9 a = mod(bitshift(a, 1), m);

10 b = bitshift(b, -1);

11 end

12 c = mod(c, m);

13end

Lines 7 and 9 are where this algorithm can overflow for 64-bit integers. Line
7 adds $c+a$ without overflow guards. This is fixed by _ModAdd_ :
$c=ModAdd(c,a,m)$. Line 9 utilizes a _bitshift_ of $a$. Bitshift by one shifts
all bits to the left by one, effectively doubling a value. For values ranging
from $2^{63}$ to $2^{64}$, a naive bitshift loses the highest bit and the
result _decreases_ rather than doubles. _ModAdd_ can resolve this overflow as
well, $a=ModAdd(a,a,m)$.

These two updates fix _ModMulitplyClassic_ ’s overflow issues. This
implementation now works correctly for all 64-bit numbers by avoiding
saturating overflows. However, this revised implementation is not yet
optimized for performance due to its large _while_ loop and absence of
vectorized MATLAB instructions.

_Vectorization_ in MATLAB is the operation on multiple elements of an array in
one calculation [7]. For example,

        >> [1, 2, 3] + 1
        ans =

             2     3     4

The scalar $1$ implicitly expands into an array then is added element-wise to
$[1,2,3]$. Vectorization can greatly improve performance due of the memory
layout of MATLAB variables: the _data_ portion of variables in MATLAB occupy a
contiguous memory buffer [6]. Vectorized operations operate on this entire
chunk of memory at once. Contiguous memory is optimized for prefetching, data
locality, and branchless execution leading to faster performance [1].

_ModMultiplyClassic_ can be improved by vectorizing its independent
calculations. Note that every $b$ value in _ModMultiplyClassic_ is independent
of all other calculations: $b$ is composed of bitshifts followed by checking
individual bits. These operations can be lifted out of the _while_ loop and
vectorized into:

    shift = bitshift(b, 0:-1:-63);
    oddIdx = find(mod(shift, 2) == 1);

All 64 bitshifts are done in one statement then vectorized division by two is
performed for this array of results to determine which elements are odd. By
removing $b$ from the _while_ loop, the loop becomes simpler, faster, and
easier to optimize further.

The final version of _ModMultiply_ contains one additional vectorized
instruction: accumulating the sequence of all intermediate $a$ values created
in the loop. By storing the entire sequence into an array of length 64, each
element can be added using MATLAB’s vectorized _sum_ function:
$c=sum(sequenceA)$. The performance gain for this technique regularly
outweighs the cost to allocate 512 bytes (64 elements * 8 bytes = 512) of
memory.

###  2.4 ModExp

_ModExp_ (Modular Exponentiation) calculates $c=a^{b}\pmod{m}$ without
overflow. The bulk of its work is handled by the two building blocks _ModAdd_
and _ModMultiply_. _ModExp_ uses exponentiation by squaring to calculate the
result in logarithmic time [11].

###  2.5 Summarizing the Miller–Rabin primality test

_ModAdd_ , _ModMultiply_ , and _ModExp_ combine to implement the Miller–Rabin
test without arbitrary-precision arithmetic. Miller–Rabin is a modern
algorithm, utilizing powerful prime number congruence relations with a careful
selection of bases to determine primality deterministically. Using these
tools, _isprime_fast_ decidedly beats the speed of _isprime_ for large
integers: for inputs in the 249 to 264 range, _isprime_fast_ ranges from 5x to
2950x faster than _isprime_ (results summarized in Fig. 2).

##  3 Improving the performance for medium scalar integers (218 to 249)

MATLAB’s _isprime_ performs well in this range, however a few helpful
techniques make it run 3-8x faster for primes and 7-30x faster for composites.

Profiling _isprime(N)_ reveals that most of its execution time is spent
sieving for primes $\leq$ _sqrt(N)_. To improve _isprime_ , _isprime_fast_
eliminates the sieve of Eratosthenes completely in favor of a simpler wheel
sieve. A wheel sieve is an increasing sequence of numbers with select small
primes removed. For example, a 2-wheel is $[1,3,5,7,9,...]$, and a 2/3-wheel
is $[1,5,7,11,13,...]$. After creating this sequence, vectorized trial
division is run against the input scalar $N$. If none of the divisions divides
evenly, $N$ must be prime. Efficiently generating a wheel sieve requires
understanding two properties of the wheel:

  1. 1.

The cyclical _differences_ between elements of the wheel

  2. 2.

The _cycle length_

For example, a 2/3-wheel begins with $[1,5,7,11,13,...]$. The difference
between consecutive elements repeats: $[4,2,4,2,...]$. Its _cycle length_
equals the product of the removed primes: $2*3=6$. Thus, a 2/3-wheel has a
difference cycle of $[4,2]$ which repeats every 6 numbers. Similarly, a
2/3/5-wheel has a difference cycle of $[6,4,2,4,2,4,6,2]$ which repeats every
$2*3*5=30$ numbers. And the 2/3/5/7-wheel used by _isprime_fast_ has a
repeating difference sequence of:

$[10,2,4,2,4,6,2,6,4,2,4,6,6,2,6,4,2,6,4,6,8,4,2,4,2,4,8,6,4,6,2,4,6,2,6,6,4,2,4,6,2,6,4,2,4,2,10,2]$
(5)

which repeats every $2*3*5*7=210$ numbers.

By knowing the sequence of differences and the cycle length, _isprime_fast_
can efficiently replicate the repeating sequence up to the length of
_sqrt(N)_. Then, to get the integers needed for trial division, the elements
are _cumulatively summed_ (_cumsum_). This is when the total sum accumulates
in each successive element. For example, to obtain the first six elements
needed for division in a 2/3/5/7-wheel, compute:

$1+cumsum([10,2,4,2,4,6])=1+[10,12,16,18,22,28]=[11,13,17,19,23,29]$ (6)

Cumulative sum of a repeating sequence is several times faster than _isprime_
’s sieve for two primary reasons. First, MATLAB excels at creating sequences
using its built-in vectorized functions such as _cumsum_ and _repmat_ (_rep_
eat _mat_ rix). Second, creating the sequence requires only one _backtrack_
(i.e., visiting a previously visited element). By comparison, the sieve of
Eratosthenes requires extra backtracks to cross-off all composite numbers.
This leads to its time complexity of _O(N log log N)_ [12] with comparatively
high constant factors compared to a simpler wheel. Creating a wheel in MATLAB
is a _linear O(N)_ operation which requires only two visits of each element:

  1. 1.

The cycle sequence is repeated up to length _sqrt(N)_

  2. 2.

These elements are added in a cumulative sum of length _sqrt(N)_

_isprime_fast_ removes the bottleneck prime number sieve and focuses on
efficient, vectorized operations. Its speedup factor for medium size scalar
inputs is typically 3x to 30x faster than _isprime_ (results summarized in
Table 4 and Fig. 2).

##  4 Improving the performance for small scalar integers (0 to 218)

For small prime scalars, _isprime_fast_ uses a 2-wheel (i.e., odd numbers)
instead of a more elaborate 2/3/5/7-wheel. Inputs in this range are so small
that the cost of repeating then summing the sequence of differences becomes
significant. Odd numbers have a constant difference between elements.
Therefore, creating a sequence of odd numbers requires only a single visit of
each element. By checking divisibility against odd numbers using integer
division, _isprime_fast_ typically beats _isprime_ by a factor of 5 to 10 for
scalar inputs in this range. Section 6.1 further emphasizes the importance of
integer division.

##  5 Symbolic Math Toolbox

The Symbolic Math Toolbox is a supplemental MATLAB product which includes a
suite of mathematical functions for equation solving, factoring, prime
numbers, and more. The toolbox supports arbitrary-precision arithmetic and is
not included with MATLAB. Both _isprime_ and _isprime_fast_ have no
dependencies on it or any other arbitrary-precision calculations.

For determining prime, the Symbolic Math Toolbox includes a function called
_sym/isprime_ which utilizes a third-party C library and yields excellent
performance for large integers (greater than about 250). For small integers
and long arrays, however, most of the execution time is spent converting data
recognizable by one programming interface into the _same_ data recognizable by
a _different_ programming interface.

_isprime_fast_ is intended to replace MATLAB’s _isprime_ and _not_
_sym/isprime_. For completeness, the performance of all three (_isprime_fast_
, _isprime_ , and _sym/isprime_) is compared.

##  6 Improving _isprime_ ’s performance for non-scalar inputs

Until now, only _scalar_ (1x1) inputs have been considered. MATLAB’s _isprime_
function also supports input arrays, matrices, and higher dimensional inputs.
Recall that MATLAB’s _isprime(N)_ implementation consists of two parts:

  1. 1.

Sieve. The sieve returns all prime numbers up to the square root of the
maximum input of the array.

  2. 2.

Trial division of primes. _isprime_ divides each element of input _N_ by all
primes which are less than that element.

These two steps have a trade-off: sieving is kept efficient by calculating
primes up to $sqrt(max(N))$ instead of $max(N)$. However, the second step
requires many divisions of those primes to determine primality. Consider these
two steps for _isprime(1:10)_ :

  1. 1.

Sieve. _floor(sqrt(10))_ is 3, so the primes to be checked against for trial
division are those $\leq 3$, which are $[2,3]$.

  2. 2.

Trial division. Dividing element-wise from 1 to 10,

    1. (a)

1 is greater than neither 2 nor 3, so there are 0 divisions

    2. (b)

2 is greater than neither 2 nor 3, so there are 0 divisions

    3. (c)

3 is greater than 2 but not 3, so there is 1 division (3 divided by 2)

    4. (d)

The elements $[4,5,6,7,8,9,10]$ are greater than both 2 and 3, so each of
these 7 elements is divided by both 2 and 3. There are 2 divisions * 7
elements = 14 total divisions

Thus, _isprime(1:10)_ performs $0+0+1+14=15$ total divisions. The number of
divisions required for incrementing sequences $[1,2,3,...,\emph{N}]$ can be
calculated by following MATLAB source. The results are in Table 1.

    p = primes(floor(sqrt(N)));
    numDivisions = sum(arrayfun(@(n) nnz(p < n), 1:N))

Input (1:N) | Number of divisions | Number of divisions / N
---|---|---
1:10 | 15 | 1.5
1:100 | 383 | 3.8
1:1000 | 10,840 | 10.8
1:10000 | 248,940 | 24.8
1:100000 | 6,490,794 | 64.9
1:1000000 | 167,923,873 | 167.9
1:10000000 | 4,459,357,131 | 445.9
1:100000000 | 122,894,263,604 | 1228.9
Table 1: Number of divisions that _isprime_ requires for incrementing
sequences.

###  6.1 Why the total number of divisions is important

Floating-point division is a comparatively slow mathematical operation on
CPUs, often 10 times slower than multiplication [3]. Demonstrated by Table 1,
the total number of divisions _accelerates_ with longer incrementing
sequences. Therefore, _isprime_ ’s performance degrades at an increasing rate
for these larger data sets, leaving enormous room for improvement.

Instrumenting _isprime_ for the array of 1 to 1,000,000 indicates that
_isprime_ spends 99.9% of its duration performing remainder division. To
improve this, _isprime_fast_ uses three techniques:

  1. 1.

Improving the speed of division

  2. 2.

Shrinking the input array

  3. 3.

Minimizing the number of divisions

####  6.1.1 Improving the speed of division

The default and most common data type in MATLAB is double [6], so _isprime_
performs floating-point division for most inputs. Floating-point division
preserves extra decimal precision which is not needed given that primality
only considers whole numbers. A simple update to use integer division
significantly improves _isprime_ ’s performance:

Calculation | Duration (sec)
---|---
isprime(1:1000000) | 8.702
isprimeIntDiv(1:1000000) | 1.377
Table 2: Comparison of floating-point division vs. integer division.

This 6.3x speedup requires only a simple data cast to _int32_. _isprimeIntDiv_
is otherwise identical to _isprime_.

####  6.1.2 Shrinking the input array

One useful optimization is to remove multiples of small primes from the input
array before additional processing. For sequences which increment by one, this
significantly reduces the size of the input array. Removing elements divisible
by 2 shrinks the array size by 50%, removing 2/3 shrinks by 66%, 2/3/5 shrinks
by 73%, and 2/3/5/7 shrinks by 77% [14, 15]. _isprime_fast_ trims $[2,3,5,7]$
which works well based on heuristics given that removing each additional prime
has diminishing returns.

####  6.1.3 Minimizing the number of divisions

Instead of division, a strategy which works well for many long arrays is to
calculate _all_ primes up to _N_ rather than the _sqrt(N)_. Having all primes
avoids division completely: _isprime_fast_ simply needs to check for existence
of each element of the input array in the set of all primes. MATLAB’s prime
number sieve returns a _sorted_ list of primes, so _binary search_ efficiently
finds which inputs are prime. MATLAB’s versatile _ismember_ function is suited
for this.

Table 3 compares the duration of incrementing sequences for _isprime_ ,
_isprimeIntDiv_ , and the final _isprime_fast_ which includes all three
optimizations.

Input |  _isprime_ (sec) |  _isprimeIntDiv_ (sec) |  _isprime_fast_ (sec) |  _isprime_fast_ speedup
---|---|---|---|---
1:10000 | 0.0185 | 0.0101 | 0.000829 | 22.3x
1:100000 | 0.3352 | 0.1174 | 0.005502 | 60.9x
1:1000000 | 8.70 | 1.794 | 0.05697 | 152.7x
1:10000000 | 139.50 | 20.86 | 0.6294 | 221.6x
1:100000000 | 2898.13 | 820.86 | 6.728 | 430.7x
Table 3: Performance of checking primality of incrementing sequences.

The performance for _isprime_fast_ scales in an approximately linearithmic
manner, _O(N lg N)_ , due to its binary search (_lg N_) of each of the _N_
elements. _isprime_ slows much more rapidly with larger input sizes due to the
accelerating number of divisions that it requires (see Table 1).

However, incrementing sequences represent only a tiny sample of all possible
input arrays. The open question is: when to choose the _sqrt_ path and when to
choose the _binary search_ path?

###  6.2 Reconciling the _sqrt(N)_ path with the _binary search_ path

There are two possible execution flows for arrays in _isprime_fast_ to be
compared:

  1. 1.

Sieving primes up to _sqrt(N)_ then integer division of primes.

  2. 2.

Sieving primes up to _N_ then binary searching for primes.

The decision of which path to take is determined through heuristics. Data
tested includes random 32-bit primes, random numbers, random odd numbers,
incrementing sequences, odd numbers in a restricted range, odd numbers with
small primes removed, pseudoprimes, random tiny and small numbers, data with a
normal distribution, and more. Based on this assortment of data, the two
biggest influences on performance are (1) the number of elements in the array,
and (2) the maximum value within the array. To determine the best fit, arrays
were constructed which have roughly _equal_ duration for each of the two
paths.

The raw performance data is noisy. Due to this volatility, _isprime_fast_ uses
a linear best fit to avoid overfitting. With a first-degree polynomial fit,
the raw data’s R2 value is low, $0.694/1$. Additionally, the _binary search_
path usage is niche: it performs superlatively for long arrays with a
relatively small maximum value but degrades quickly if misused on shorter
arrays with a high maximum value. Conversely, the _sqrt_ path (which is what
MATLAB’s _isprime_ always uses) performs adequately for all varieties of
inputs. Therefore, from the raw data, only the data points with comparatively
low maximum values were selected for the final curve fitting. This errors on
the side of caution: binary search is used only when the input data metrics
abundantly work in its favor.

When considering this subset of data points, the data forms two distinct
clusters: one for a small number of elements (less than about 30,000) and one
for a larger number of elements. Figure 1 shows the inflection points where
both the _sqrt_ and _binary search_ paths have equal duration.

Figure 1: Arrays of differing length and maximum values when the _sqrt_ path
and _binary search_ path have equal duration
$(A)$ Low number of elements, $R^{2}=0.99275$. $(B)$ High number of elements,
$R^{2}=0.91269$.

For a small number of elements in the input array, the linear best fit is
approximately:

$275*numElements-100000$ (7)

For the longer data, the linear fit is:

$613*numElements-200000000$ (8)

When these equations return a value _smaller_ than the maximum input value,
the _sqrt_ path is taken. Otherwise, binary search is selected. These
equations perform well for many categories of input arrays, are backed by
performance tests, and can readily be adjusted based on new results.

##  7 How the results were obtained

###  7.1 Software and hardware

All results were obtained using MATLAB R2020b on Windows 10, Intel® Core™
i7-9750H CPU @ 2.60GHz, 16 GB RAM. Results were also verified on:

  * •

Windows 10, Intel® Xeon® W-2133 CPU @ 3.60GHz, 64 GB RAM

  * •

Debian GNU/Linux 10 (buster), Intel® Xeon® Gold 6140 CPU @ 2.30GHz, 51 GiB RAM

  * •

macOS Big Sur, Quad-Core Intel® Core i3 @ 3.60GHz, 32 GB RAM

###  7.2 Reproducibility

The results on both Windows hardware were mostly consistent within a few
percentage points. Linux and Mac fluctuated more but _isprime_fast_ continued
to outperform _isprime_ strongly for the suite of performance tests in the
GitHub repository.

###  7.3 Source Code

The _isprime_fast_ MATLAB source is available on GitHub,
https://github.com/tnear/isprime_fast.

##  8 Results

###  8.1 Scalar performance comparison

Table 4 summarizes the performance of _isprime_fast_ , MATLAB’s _isprime_ ,
and Symbolic’s _sym/isprime_ for prime scalars. The first column, bit-size, is
the number of bits of the numeric input. For _isprime_ and _isprime_fast_ ,
double-precision is used for inputs $\leq$ 253 and unsigned 64-bit integer for
inputs > 253 (IEEE 754 double-precision cannot represent odd numbers higher
than 253). The symbolic _sym_ data type is used for _sym/isprime_. The timings
listed are the average duration of _100 random prime_ number inputs with the
specified number of bits.

Bit-size |  _isprime_fast_ (sec) |  _isprime_ (sec) |  | _isprime_fast_
---
speedup
_sym/isprime_ (sec) |  | _isprime_fast_
---
speedup
4 | 7.7654e-07 | 4.7186e-06 | 6.07x | 0.002166 | 2789x
8 | 1.0120e-06 | 5.7043e-06 | 5.63x | 0.00209 | 2066x
16 | 1.8106e-06 | 1.1762e-05 | 6.49x | 0.002196 | 1213x
24 | 5.8244e-06 | 4.3816e-05 | 7.52x | 0.008325 | 1429x
32 | 5.3932e-05 | 4.1277e-04 | 7.65x | 0.0124 | 230x
36 | 0.0002828 | 0.0011006 | 3.89x | 0.012055 | 42.6x
40 | 0.001766 | 0.004343 | 2.45x | 0.016207 | 9.17x
44 | 0.005818 | 0.01744 | 3.00x | 0.015349 | 2.64x
48 | 0.02718 | 0.07891 | 2.90x | 0.01777 | 0.653x
50 | 0.0326 | 0.1901 | 5.83x | 0.0211 | 0.648x
52 | 0.0431 | 0.3973 | 9.21x | 0.0244 | 0.566x
56 | 0.0552 | 1.8212 | 32.9x | 0.0226 | 0.410x
60 | 0.0755 | 7.4418 | 98.6x | 0.0280 | 0.370x
64 | 0.1194 | 31.8824 | 266x | 0.0287 | 0.234x
Table 4: Performance comparison of _prime_ scalar inputs of various bit-sizes.

Table 4 compares the performance of prime inputs only. Figure 2 summarizes
_isprime_fast_ ’s performance improvement over _isprime_ for prime numbers,
odd numbers, and random numbers.

Figure 2: _isprime_fast_ speedup factor compared to _isprime_ for various
scalar input bit-sizes.

The solid line on the bottom denotes random prime numbers, the dashed line is
random odd numbers (complete data in Table 6), and the dash-dot line at the
top is random numbers (complete data in Table 7).

The improvement of _isprime_fast_ over _isprime_ ranges from about 2.4x to
260x faster for random prime numbers depending on the bit-size of the input.
For random odd numbers, _isprime_fast_ executes about 5.4x to 1600x faster
than _isprime_. And testing random numbers widens the gap further:
_isprime_fast_ ’s improvement ranges from 7.3x to 2950x faster.

The comparison between _isprime_ and _sym/isprime_ in Table 4 illustrates the
importance of data types. _isprime_fast_ is thousands of times faster for
small integers: symbolic instructions leave the MATLAB language and incur a
noticeable cost for computations which otherwise take microseconds. However,
for large primes, _sym/isprime_ is 2x to 4x faster than _isprime_fast_ and 9x
to 1000x faster than _isprime_. The time that _sym/isprime_ takes to convert
symbolic MATLAB data into a format recognizable to the C library is more than
compensated by that library’s faster calculation.

###  8.2 Array performance comparison

Table 5 summarizes the performance comparison results for _isprime_fast_ ,
_isprime_ , and _sym/isprime_ for a variety of array input data types, array
lengths, bit-sizes, and data distributions.

Input array |  | _isprime_fast_
---
(sec)
| _isprime_
---
(sec)
| _isprime_fast_
---
speedup
| _sym/isprime_
---
(sec)
| _isprime_fast_
---
speedup
1:10000000 (10M) | 0.6294 | 139.50 | 221x | 79325* | 126032x
| 1609477200:1:(
---
1609477200 +
3600*24*365)
(Unix seconds
in the 2021
calendar year)
12.26 | 3286.55 | 268x | 26185 | 2136x
| 200 random
---
int16 primes
5.0952e-05 | 2.2113e-04 | 4.34x | 0.0245 | 480x
| 1000x1000 matrix
---
of random 32-bit
integers
3.927 | 16.647 | 4.24x | 529.9 | 134x
| Normal distribution
---
of 1 million primes
with mean 216
0.1623 | 5.370 | 33.0x | 7614.4 | 46915x
| Normal distribution
---
of 100,000 random
odd numbers with
mean 232
2.689 | 19.603 | 7.28x | 144.6 | 53.7x
| Normal distribution
---
of 1,000 random
numbers with
mean 253
5.07 | 74.32 | 14.6x | 0.3370 | 0.0665x
| The largest 52-bit
---
prime repeated
1,000 times
58.270 | 172.516 | 2.96x | 18.445 | 0.3165x
| The 100 largest
---
pseudoprimes >253
with factors >20M
0.4852 | 7.4524 | 15.3x | 0.01160 | 0.0239x
| The 100 largest
---
64-bit primes
0.1442 | 206.56 | 1432x | 0.0284 | 0.197x
| The 5000 largest
---
64-bit odd integers
76.201 | 8816.8 | 115.7x | 7.951 | 0.1043x
Table 5: Performance comparison of non-scalar inputs for _isprime_fast_
_*estimated_

Based on Table 5, _isprime_fast_ performs strongly over _isprime_ for arrays
containing incrementing sequences, prime/odd/composite numbers, normal
distributions, random small/medium/large primes, pseudoprimes, and
short/medium/long arrays. The performance difference is greatest for large
primes and long arrays where _isprime_fast_ often operates hundreds of times
faster. The difference is smallest for small integer data types where
_isprime_fast_ is only about 2-4x faster. This is due to _isprime_ executing
faster integer (instead of floating-point) division— _isprime_fast_ has fewer
opportunities for optimizations over it.

The metrics for Symbolic’s _sym/isprime_ reflect the explanation after Table 4
about symbolic data. The cost per element to convert to a format recognizable
by the external library is massive for long arrays: two of the performance
benchmarks indicate that _isprime_fast_ is over 46,000x faster than
_sym/isprime_. But for large primes inside short arrays, _sym/isprime_ handily
beats both _isprime_ and _isprime_fast_.

##  9 Discussion

The largest improvement of _isprime_fast_ over MATLAB’s _isprime_ is for
64-bit integers where _isprime_fast_ runs 250x to 3000x faster. The modern
Miller—Rabin algorithm for 64-bit numbers is tremendously fast compared to
_isprime_ ’s classic prime number sieve followed by prime factor division. One
reason for this great contrast is because _isprime_ ’s implementation predates
MATLAB’s support for 64-bit integers: the _int64_ and _uint64_ data types were
not supported until 2004 [8]. Before then, MATLAB only supported the double-
precision data type which has a maximum odd value of $2^{53}-1$. This is
significant because _isprime_ completes in one second for 53-bit numbers but
its duration balloons to over 30 seconds for 64-bit integers.

The largest area for improvement of _isprime_fast_ is for scalar inputs in the
40- to 48-bit range. Figure 2 shows a dip in performance for input scalars in
this range: _isprime_fast_ only executes two to three times faster than
_isprime_. _isprime_fast_ uses a 2/3/5/7-wheel sieve for inputs in this range.
For 48-bit integers, this wheel is over 3 million elements long and therefore
requires equally as many divisions to determine primality. An alternative
algorithm likely could remedy this performance depression.

##  10 Conclusions

MATLAB’s _isprime_ implementation is stable, effective, and depended on by
many mathematics enthusiasts and researchers. _isprime_ provides fair
performance for small to medium size scalars and short arrays, but less than
optimal performance otherwise. _isprime_fast_ , a proposed replacement for
_isprime_ , improves upon _isprime_ for all 264 scalar inputs, ranging from
2.4x to 2950x faster. _isprime_fast_ also surpasses _isprime_ for all arrays
which were tested, normally ranging from 2x to 1400x faster. _isprime_fast_
efficiently creates wheel sieves, implements the Miller—Rabin test using
modular arithmetic strategies, minimizes divisions, and heavily relies on
vectorized operations. _isprime_fast_ employs the power of the MATLAB
language: it requires no arbitrary-precision arithmetic, C/C++ source code, or
external libraries. Using the techniques detailed here, MATLAB’s esteemed
_isprime_ function can be improved considerably.

###  10.1 Future research

_isprime_fast_ still has ample for improvement and additional research. The
findings here are just the beginning into the large field of prime number
research. _isprime_fast_ uses deterministic Miller–Rabin for large primes.
While this is hundreds of times faster than _isprime_ ’s sieve with trial
division for 64-bit integers, it is not the fastest prime number checker. The
Baillie–PSW primality test, for example, utilizes one iteration of
Miller–Rabin with base 2 followed by a strong Lucas probable prime test [9].
Implementing the Lucas test without arbitrary-precision arithmetic requires
considerable reliance on _ModExp_ , _ModMultiply_ , and _ModAdd_ plus
additional modular arithmetic considerations to avoid overflowing 64-bit
integers. The results are promising: implemented in MATLAB, Baillie–PSW often
runs 1.5x to 2x faster than Miller–Rabin. However, due to the fluctuating
nature of Jacobi symbol distributions and the complexity of the Lucas test,
the MATLAB implementation still needs additional testing.

## References

  * [1] U. Drepper.  What every programmer should know about memory.  Red Hat, Inc., pages 12,14,56, 2007.
  * [2] GeeksforGeeks.  How to avoid overflow in modular multiplication?  https://www.geeksforgeeks.org/how-to-avoid-overflow-in-modular-multiplication/, 2021\.
  * [3] I. Hare.  Not all cpu operations are created equal. infographics: Operation costs in cpu clock cycles.  http://ithare.com/infographics-operation-costs-in-cpu-clock-cycles/, 2016\.
  * [4] B. Lynn.  Number theory – primality tests.  https://crypto.stanford.edu/pbc/notes/numbertheory/millerrabin.html.
  * [5] MathWorks.  Determine which array elements are prime.  https://www.mathworks.com/help/matlab/ref/isprime.html, 2006.
  * [6] MathWorks.  Matlab data.  https://www.mathworks.com/help/matlab/matlab_external/matlab-data.html, 2021\.
  * [7] MathWorks.  Using vectorization.  https://www.mathworks.com/help/matlab/matlab_prog/vectorization.html, 2021\.
  * [8] C. Moler and J. Little.  A history of matlab.  Proc. ACM Program. Lang. 4, HOPL, Article 81, page 24, 2020.
  * [9] C. Pomerance, J. L. Selfridge, and S. S. Wagstaff.  The pseudoprimes to 25*109.  Mathematics of Computation, 35:1003–1026, 1980.
  * [10] M. Rabin.  Probabilistic algorithm for testing primality.  Journal of Number Theory, pages 128–138, 1980.
  * [11] B. Schneier.  Applied Cryptography: Protocols, Algorithms, and Source Code in C.  Wiley, 1996.
  * [12] R. Sedgewick.  Algorithms in C++.  Addison-Wesley, 1992.
  * [13] J. Sinclair.  Deterministic variants of the miller-rabin primality test.  https://miller-rabin.appspot.com/, 2011.
  * [14] N. J. A. Sloane.  a(0) = 1; for n > 0, a(n) = (prime(n)-1)*a(n-1).  https://oeis.org/A005867.
  * [15] N. J. A. Sloane.  Primorial numbers: product of first n primes.  https://oeis.org/A002110.

##  11 Appendix tables

Bit-size |  _isprime_fast_ (sec) |  _isprime_ (sec) |  | _isprime_fast_
---
speedup
_sym/isprime_ (sec) |  | _isprime_fast_
---
speedup
4 | 7.9037e-07 | 4.5998e-06 | 5.82x | 0.00211 | 2668x
8 | 9.8097e-07 | 5.4535e-06 | 5.56x | 0.0020 | 2084x
16 | 1.4468e-06 | 9.6701e-06 | 6.68x | 0.0020 | 1409x
24 | 2.9905e-06 | 3.9433e-05 | 13.1x | 0.00299 | 1000x
32 | 1.9783e-05 | 3.5878e-04 | 18.1x | 0.00298 | 150x
36 | 2.3884e-05 | 4.0707e-04 | 17.0x | 0.002895 | 121x
40 | 7.2903e-04 | 0.0040 | 5.48x | 0.0027 | 3.49x
44 | 0.002690 | 0.01878 | 6.98x | 0.003328 | 1.23x
48 | 0.01062 | 0.07565 | 7.1x | 0.002771 | 0.260x
50 | 0.00636 | 0.1616 | 25x | 0.00279 | 0.439x
52 | 0.007706 | 0.39378 | 51x | 0.002832 | 0.367x
56 | 0.0095 | 1.8229 | 192x | 0.0044 | 0.463x
60 | 0.01316 | 7.4327 | 564x | 0.003645 | 0.277x
64 | 0.0196 | 31.7576 | 1620x | 0.0034 | 0.174x
Table 6: Performance comparison of random _odd_ scalar inputs of various bit-sizes. Bit-size |  _isprime_fast_ (sec) |  _isprime_ (sec) |  | _isprime_fast_
---
speedup
_sym/isprime_ (sec) |  | _isprime_fast_
---
speedup
4 | 5.4788e-07 | 4.0505e-06 | 7.39x | 0.002046 | 3734x
8 | 7.1275e-07 | 5.4507e-06 | 7.64x | 0.002103 | 2950x
16 | 9.0833e-07 | 9.7864e-06 | 10.7x | 0.002008 | 2210x
24 | 1.8587e-06 | 4.0482e-05 | 21.7x | 0.002289 | 1231x
32 | 1.1610e-05 | 3.7037e-04 | 31.9x | 0.002584 | 222x
36 | 6.9472e-05 | 0.001146 | 16.5x | 0.002507 | 36.0x
40 | 3.3226e-04 | 0.004121 | 12.4x | 0.002913 | 8.76x
44 | 0.001256 | 0.01854 | 14.7x | 0.002359 | 1.87x
48 | 0.004622 | 0.06985 | 15.1x | 0.002308 | 0.4993x
50 | 0.003334 | 0.1755 | 52.6x | 0.002601 | 0.7801x
52 | 0.003863 | 0.3647 | 94.4x | 0.002498 | 0.6467x
56 | 0.005152 | 1.7537 | 340x | 0.003188 | 0.6189x
60 | 0.006451 | 7.2382 | 1122x | 0.002881 | 0.4466x
64 | 0.01126 | 33.2620 | 2953x | 0.002531 | 0.2247x
Table 7: Performance comparison of _random_ scalar inputs of various bit-
sizes.


simplify comparing forces we will represent $(\frac{dE}{dx})_{f}$ as
$m^{2}V(1+z)^{5/2}\mathcal{F}$. $\mathcal{F}$ represents the mostly $z$, $m$
and $V$ independent component of $(\frac{dE}{dx})_{f}$

If $\mathcal{F}z^{5/2}M>H(z)$ at the regular decoupling redshift $z_{d}$, then
the EMBHs will continue moving apart with the Hubble flow until
$\frac{2G}{x^{2}(1+z)^{1/2}\mathcal{F}}=\frac{H(z)x}{(1+z)}$. During radiation
domination this gives a decoupling redshift $z_{b}$ of

$(1+z_{b})=\Big{(}\frac{2G}{H_{0}\sqrt{\Omega_{r}}\mathcal{F}}\Big{)}^{2/3}x^{-2}=k_{b}x^{-2}$
(91)

Where we will denote
$k_{b}\equiv\Big{(}\frac{H_{0}\sqrt{\Gamma_{r}}\mathcal{F}}{2G}\Big{)}^{2/3}$,
and where $\Omega_{r}$ is the fraction of the critical energy density composed
of radiation at $z=0$ and $z_{f}$ denotes the redshift at which frictional
forces become small compared to Hubble:

$(1+z_{f})=\frac{\mathcal{F}M^{2}}{H_{0}\sqrt{\Omega_{r}}}$ (92)

Binaries must decouple by matter radiation equality because the ratio between
the free fall time and the Hubble time becomes constant during matter
domination, so those binaries that could not decouple by matter radiation
equality will not gain the ability to do so during matter domination,

The semimajor axis of the binary $a_{0}$ is just $\frac{x}{2(1+z_{d})}$. The
eccentricity $e_{0}$ of the binary is set by interactions with nearest EMBH
outside of the binary. One of the EMBHs in the binary gets pulled out of it’s
free fall path by the next nearest EMBH, while the other does not strongly
feel any forces from the next nearest charge. The distance the in-falling EMBH
gets pulled out of its free fall path will set the semi-minor axis for the
binary system $b$. This can be calculated by multiplying the force exerted by
the EMBH nearest the binary with the square of the time it takes the perturbed
EMBH to fall a distance $a_{0}$. Setting the comoving distance to the nearest
EMBH outside of the binary to $y$ this approach gives a $b$ of

$2b=\frac{1}{2}\Big{(}\frac{2GM(1+z_{d})^{2}}{y^{2}}\Big{)}\Big{(}\frac{x^{3}}{2GM(1+z_{d})^{3}}\Big{)}=\frac{x^{2}}{y^{2}}a$
(93)

This gives a starting eccentricity of

$e_{0}=\sqrt{1-\frac{x^{4}}{4y^{4}}}$ (94)

The probability that two oppositely charged EMBHs are separated by a co-moving
distance between $x$ and $x+dx$ and that the next nearest EMBH of any sign
charge is between $y$ and $y+dy$.

$dP=\frac{1}{2}(4\pi n)^{2}x^{2}y^{2}dxdy$ (95)

where $n$ is the comoving number density of EMBH. The factor of $\frac{1}{2}$
arises because each EMBH can only pair with an oppositely charged partner, so
the average distance between "pairable" EMBHs is $2^{1/3}$ times the average
distance between all EMBH. $x$ is by definition less than $y$, because a
binary would just form with the EMBHs separated by $y$ otherwise. $y$ can be
at most the average separation between EMBH. This defines $y_{max}$

$y_{max}=\frac{3M}{4\pi\rho_{dm}f}$ (96)

This can be re-written using (89) and (94) to give a probability distribution
for the orbital parameters of all of the binaries. That can then be related to
the probability that a binary merges at a specific time $t_{m}$ by (88) which
is derived in Section B.1. Rearranging (88) gives $a_{0}(t_{m},e_{0})$, the
starting semi-major axis that corresponds to a binary with eccentricity
$e_{0}$ that would merge in $t_{m}$.

$a_{0}=\frac{M^{2/3}t_{m}^{1/3}}{\left(\frac{0.034}{8G^{2}}\right)^{1/3}(1-e_{0}^{2})}$
(97)

All binaries that merge at $t_{m}$ have parameters that lie along the curve in
$e-a$ space $a_{0}(t_{m},e_{0})$. With a bit of algebra, one finds that in the
limit where $e\sim 1$ the probability that an EMBH binary merges between
$t_{m}$ and $t_{m}+dt$ is

$dP=\frac{(4\pi
n)^{2}}{24}\frac{M^{5/2}k_{d}^{3/2}}{\left(\frac{0.034}{8G^{2}}t_{m}\right)^{1/2}}\frac{e_{0}}{(1-e_{0}^{2})^{13/4}}dtde$
(98)

The merger rate per EMBH, $\Gamma_{m}$ at a specific time $t_{m}$ can be found
then by integrating $\frac{dP}{dt}$ over all eccentricities a binary merging
at $t_{m}$ could have.

$\Gamma_{m}=\frac{dP}{dt}=\frac{(4\pi\rho_{DM}f)^{2}}{27*4}\frac{M^{1/2}k_{d}^{3/2}}{\left(\frac{0.034}{8G^{2}}t_{m}\right)^{1/2}}\Bigg{(}\frac{1}{(1-e_{max}^{2})^{9/4}}-\frac{1}{(1-e_{min}^{2})^{9/4}}\Bigg{)}$
(99)

Multiplying $\Gamma_{m}$ by $n_{bh}$, the proper number density of EMBHs at
$t_{m}$, would give the rate at which EMBH mergers happen at $t_{m}$. The
bounds used in these integrals $e_{max}$ and $e_{min}$ are the largest and
smallest eccentricities a binary merging at $t_{m}$ could have.

We find the limits on eccentricities as follows. We define $e_{upper}$ as the
highest eccentricity a binary with semi-major axis $a$ can have. This limit
exists because $y$ cannot be greater than the average spacing between EMBHs
$y_{max}=\Big{(}\frac{3M}{4\pi\rho_{dm}f}\Big{)}^{1/3}$. This gives

$e_{upper}=\sqrt{1-\frac{a}{4k_{d}y_{max}^{4}}}$ (100)

Then, $e_{max}$ is set by where the curve $a(e_{0},t_{m})$, intersects
$e_{upper}$, or $a_{eq}$, the semi-major axis of a binary that decouples at
the latest allowed decoupling redshift $z_{eq}$.

$e_{max}=\begin{cases}\sqrt{1-\left(\frac{M^{2/3}t_{m}^{1/3}}{4k_{d}\left(\frac{0.034}{8G^{2}}\right)^{1/3}y_{m}^{4}}\right)^{1/2}}&t_{m}\leq
T_{c}\\\
\sqrt{1-\left(\frac{k_{d}M^{2}t_{m}(1+z_{eq})^{4}}{\left(\frac{0.034}{8G^{2}}\right)}\right)^{1/3}}&t_{m}>T_{c}\end{cases}$
(101)

where $T_{c}$ is

$T_{c}=\frac{\frac{0.034}{8G^{2}}}{2^{6}M^{2}k_{d}^{5}y_{m}^{12}(1+z_{eq})^{8}}$
(102)

Finally, $e_{min}$ is the lowest eccentricity a binary that merges in $t_{m}$
and decouples after $z_{f}$ can have. It is 0 when all binaries that merge in
$t_{m}$ decouple after $z_{f}$.

$e_{min}=\begin{cases}\sqrt{1-\left(\frac{k_{d}M^{2}t_{m}(1+z_{f})^{4}}{\left(\frac{0.034}{8G^{2}}\right)}\right)^{1/3}}&t_{m}>T_{f}\\\
0&t_{m}<T_{f}\end{cases}$ (103)

where $T_{f}$ is

$T_{f}=\frac{\left(\frac{0.034}{8G^{2}}\right)}{M^{2}k_{d}(1+z_{f})^{4}}$
(104)

Here $e_{min}$ is determined in the limit that $e_{min}\sim 1$, this is
appropriate because the majority of binaries merging today or earlier are
highly eccentric.

## References

  * (1) G. ’t Hooft, _Magnetic monopoles in unified gauge theories_ , _Nuclear Physics B_ 79 (1974) 276.
  * (2) A.M. Polyakov, _Particle Spectrum in the Quantum Field Theory_ , _JETP Lett._ 20 (1974) 194.
  * (3) MACRO collaboration, _Final results of magnetic monopole searches with the MACRO experiment_ , _Eur. Phys. J. C_ 25 (2002) 511 [hep-ex/0207020].
  * (4) S. Orito, H. Ichinose, S. Nakamura, K. Kuwahara, T. Doke, K. Ogura et al., _Search for supermassive relics with a 2000- ${\mathrm{m}}^{2}$ array of plastic track detectors_, _Phys. Rev. Lett._ 66 (1991) 1951.
  * (5) S. Balestra, S. Cecchini, M. Cozzi, M. Errico, F. Fabbri, G. Giacomelli et al., _Magnetic monopole search at high altitude with the slim experiment_ , _The European Physical Journal C_ 55 (2008) 57–63.
  * (6) A.H. Guth, _Inflationary universe: A possible solution to the horizon and flatness problems_ , _Phys. Rev. D_ 23 (1981) 347.
  * (7) J. Preskill, _Cosmological Production of Superheavy Magnetic Monopoles_ , _Phys. Rev. Lett._ 43 (1979) 1365.
  * (8) B.J. Carr and S. Hawking, _Black holes in the early Universe_ , _Mon. Not. Roy. Astron. Soc._ 168 (1974) 399.
  * (9) A.H. Guth and S.H.H. Tye, _Phase transitions and magnetic monopole production in the very early universe_ , _Phys. Rev. Lett._ 44 (1980) 631.
  * (10) I.J. Araya, M.E. Rubio, M.S. Martín, F.A. Stasyszyn, N.D. Padilla, J. Magaña et al., _Magnetic field generation from primordial black hole distributions_ , _Monthly Notices of the Royal Astronomical Society_ (2021) .
  * (11) J. Callan, Curtis G. and E. Witten, _Monopole Catalysis of Skyrmion Decay_ , _Nucl. Phys. B_ 239 (1984) 161.
  * (12) V. Rubakov, _Adler-Bell-Jackiw Anomaly and Fermion Number Breaking in the Presence of a Magnetic Monopole_ , _Nucl. Phys. B_ 203 (1982) 311.
  * (13) Y. Bai, J. Berger, M. Korwar and N. Orlofsky, _Phenomenology of magnetic black holes with electroweak-symmetric coronas_ , 2020.
  * (14) D. Ghosh, A. Thalapillil and F. Ullah, _Astrophysical hints for magnetic black holes_ , 2020.
  * (15) J. Maldacena, _Comments on magnetic black holes_ , 2004.06084.
  * (16) P.A.M. Dirac, _Quantised Singularities in the Electromagnetic Field_ , _Proceedings of the Royal Society of London Series A_ 133 (1931) 60.
  * (17) S.W. Hawking, _Black hole explosions?_ , _Nature_ 248 (1974) 30.
  * (18) D.L.J. Ho and A. Rajantie, _The electroweak sphaleron in a strong magnetic field_ , 2020.
  * (19) E. Poisson and W. Israel, _Inner-horizon instability and mass inflation in black holes_ , _Phys. Rev. Lett._ 63 (1989) 1663.
  * (20) D.E. Kaplan and S. Rajendran, _Firewalls in general relativity_ , _Physical Review D_ 99 (2019) .
  * (21) N. Lehner, B. Wakker and B. Savage, _CII radiative cooling of the diffuse gas in the Milky Way_ , _Astrophys. J._ 615 (2004) 767 [astro-ph/0407363].
  * (22) N. Meyer-Vernet, _Energy loss by slow magnetic monopoles in a thermal plasma_ , _The Astrophysical Journal_ 290 (1985) 21.
  * (23) C.F. McKee and J.P. Ostriker, _A theory of the interstellar medium: three components regulated by supernova explosions in an inhomogeneous substrate._ , _The Astrophysical Journal_ 218 (1977) 148.
  * (24) F.X. Timmes and S.E. Woosley, _The Conductive Propagation of Nuclear Flames. I. Degenerate C + O and O + NE + MG White Dwarfs_ , _The Astrophysical Journal_ 396 (1992) 649.
  * (25) R. Janish, V. Narayan and P. Riggins, _Type Ia supernovae from dark matter core collapse_ , _Phys. Rev._ D100 (2019) 035008 [1905.00395].
  * (26) M.A. Fedderke, P.W. Graham and S. Rajendran, _White dwarf bounds on charged massive particles_ , _Physical Review D_ 101 (2020) .
  * (27) M.A. Hollands, P.-E. Tremblay, B.T. Gänsicke, M.E. Camisassa, D. Koester, A. Aungwerojwit et al., _An ultra-massive white dwarf with a mixed hydrogen–carbon atmosphere as a likely merger remnant_ , _Nature Astronomy_ 4 (2020) 663–669.
  * (28) P. Dufour, S. Blouin, S. Coutu, M. Fortin-Archambault, C. Thibeault, P. Bergeron et al., _The Montreal White Dwarf Database: A Tool for the Community_ , in _20th European White Dwarf Workshop_ , P.E. Tremblay, B. Gaensicke and T. Marsh, eds., vol. 509 of _Astronomical Society of the Pacific Conference Series_ , p. 3, Mar., 2017 [1610.00986].
  * (29) L. Siess, _Evolution of massive AGB stars. II. model properties at non-solar metallicity and the fate of Super-AGB stars_ , _Astronomy and Astrophysics_ 476 (2007) 893.
  * (30) S.P. Ahlen and K. Kinoshita, _Calculation of the stopping power of very-low-velocity magnetic monopoles_ , _Phys. Rev. D_ 26 (1982) 2347.
  * (31) A.Y. Potekhin, D.A. Baiko, P. Haensel and D.G. Yakovlev, _Transport properties of degenerate electrons in neutron star envelopes and white dwarf cores_ , _Astron. Astrophys._ 346 (1999) 345 [astro-ph/9903127].
  * (32) F. D’Antona and I. Mazzitelli, _Cooling of white dwarfs_ , _Ann. Rev. Astron. Astrophys._ 28 (1990) 139.
  * (33) S. Bagnulo and J.D. Landstreet, _Searching for the weakest detectable magnetic fields in white dwarfs. Highly-sensitive measurements from first VLT and WHT surveys_ , _Astronomy and Astrophysics_ 618 (2018) A113 [1807.09649].
  * (34) A. Kawka, _The properties and origin of magnetic fields in white dwarfs_ , _Contributions of the Astronomical Observatory Skalnate Pleso_ 48 (2018) 228 [1801.05602].
  * (35) J. Braithwaite and H.C. Spruit, _A fossil origin for the magnetic field in a stars and white dwarfs_ , _Nature_ 431 (2004) 819–821.
  * (36) J. Nordhaus, S. Wellons, D.S. Spiegel, B.D. Metzger and E.G. Blackman, _Formation of high-field magnetic white dwarfs from common envelopes_ , _Proceedings of the National Academy of Sciences_ 108 (2011) 3135–3140.
  * (37) Fermi-LAT collaboration, _The spectrum of isotropic diffuse gamma-ray emission between 100 MeV and 820 GeV_ , _Astrophys. J._ 799 (2015) 86 [1410.3696].
  * (38) T. Ukwatta, D. Stump, J. Linnemann, J. MacGibbon, S. Marinelli, T. Yapici et al., _Primordial black holes: Observational characteristics of the final evaporation_ , _Astroparticle Physics_ 80 (2016) 90–114.
  * (39) K.V. Berghaus, M.D. Diamond and D. Kaplan, _Decays of Long-Lived Relics and Their Signatures at IceCube_ , _JHEP_ 05 (2019) 145 [1811.04939].
  * (40) A.A. Zdziarski and R. Svensson, _Absorption of X-Rays and Gamma Rays at Cosmological Distances_ , _The Astrophysical Journal_ 344 (1989) 551.
  * (41) M. Lewis, K. Freese and G. Tarle, _Protogalactic extension of the Parker bound_ , _Phys. Rev. D_ 62 (2000) 025002 [astro-ph/9911095].
  * (42) F.C. Adams, M. Fatuzzo, K. Freese, G. Tarlé, R. Watkins and M.S. Turner, _Extension of the parker bound on the flux of magnetic monopoles_ , _Phys. Rev. Lett._ 70 (1993) 2511.
  * (43) M.S. Turner, E.N. Parker and T. Bogdan, _Magnetic Monopoles and the Survival of Galactic Magnetic Fields_ , _Phys. Rev. D_ 26 (1982) 1296.
  * (44) A. Bonafede, F. Vazza, M. Brüggen, M. Murgia, F. Govoni, L. Feretti et al., _Measurements and simulation of Faraday rotation across the Coma radio relic_ , _Mon. Not. Roy. Astron. Soc._ 433 (2013) 3208 [1305.7228].
  * (45) M. Kierdorf, R. Beck, M. Hoeft, U. Klein, R.J. van Weeren, W.R. Forman et al., _Relics in galaxy clusters at high radio frequencies_ , _Astronomy and Astrophysics_ 600 (2017) A18.
  * (46) J. Derkaoui, G. Giacomelli, T. Lari, G. Mandrioli, M. Ouchrif, L. Patrizii et al., _Energy losses of magnetic monopoles and dyons in scintillators, streamer tubes and nuclear track detectors_ , _Astropart. Phys._ 10 (1999) 339.
  * (47) J. Arafune and M. Fukugita, _Velocity-dependent factors for the rubakov process for slowly moving magnetic monopoles in matter_ , _Phys. Rev. Lett._ 50 (1983) 1901.
  * (48) L. Bracci and G. Fiorentini, _Interactions of Magnetic Monopoles With Nuclei and Atoms: Formation of Bound States and Phenomenological Consequences_ , _Nucl. Phys. B_ 232 (1984) 236.
  * (49) K. Olaussen, H.A. Olsen, I. Overbo and P. Osland, _PROTON CAPTURE BY MAGNETIC MONOPOLES_ , _Phys. Rev. Lett._ 52 (1984) 325.
  * (50) L.G. Althaus, A.H. Córsico and F. De Gerónimo, _Effect of Coulomb diffusion of ions on the pulsational properties of DA white dwarfs_ , _Astronomy and Astrophysics_ 644 (2020) A55 [2010.11131].
  * (51) K. Freese, _DO MONOPOLES KEEP WHITE DWARFS HOT?_ , _Astrophys. J._ 286 (1984) 216.
  * (52) K.S. Thorne and R.D. Blandford, _Modern Classical Physics: Optics, Fluids, Plasmas, Elasticity, Relativity, and Statistical Physics_ (2017).
  * (53) S.E. Woosley, S. Wunsch and M. Kuhlen, _Carbon ignition in type ia supernovae: An analytic model_ , _The Astrophysical Journal_ 607 (2004) 921–930.
  * (54) D.L. Turcotte and G. Schubert, _Geodynamics - 2nd Edition_ (2002), 10.2277/0521661862.
  * (55) F. D’Antona and I. Mazzitelli, _Cooling of white dwarfs._ , _Annual Review of Astron and Astrophys_ 28 (1990) 139.
  * (56) N.C. Hambly, S.J. Smartt and S.T. Hodgkin, _WD 0346+246: A Very Low Luminosity, Cool Degenerate in Taurus_ , _The Astrophysical Journall_ 489 (1997) L157.
  * (57) J. Liebert, C.C. Dahn and D.G. Monet, _The Luminosity Function of White Dwarfs_ , _The Astrophysical Journal_ 332 (1988) 891.
  * (58) M. Kilic, _Cool white dwarfs_ , in _White Dwarf Atmospheres and Circumstellar Environments_ , pp. 25–52, John Wiley ’|&’ Sons, Ltd (2011), DOI [https://onlinelibrary.wiley.com/doi/pdf/10.1002/9783527636570.ch2].
  * (59) M. Kilic, J.R. Thorstensen, P.M. Kowalski and J. Andrews, _11-12 gyr old white dwarfs 30 pc away_ , _Monthly Notices of the Royal Astronomical Society: Letters_ 423 (2012) L132–L136.
  * (60) E.W. Kolb and M.S. Turner, _Limits from the Soft x-Ray Background on the Temperature of Old Neutron Stars and on the Flux of Superheavy Magnetic Monopoles_ , _Astrophys. J._ 286 (1984) 702.
  * (61) M. Kutschera, _Neutron stars: formation and structure._ , _Acta Physica Polonica B_ 29 (1998) 25 [astro-ph/9801235].
  * (62) J.A. Harvey, _Magnetic monopoles in neutron stars_ , _Nuclear Physics B_ 236 (1984) 255.
  * (63) D. McCammon, D.N. Burrows, W.T. Sanders and W.L. Kraushaar, _The soft X-ray diffuse background._ , _The Astrophysical Journal_ 269 (1983) 107.
  * (64) B. Haskell and A. Sedrakian, _Superfluidity and Superconductivity in Neutron Stars_ , vol. 457, pp. 401–454 (2018), DOI [1709.10340].
  * (65) J.A. Harvey, M.A. Ruderman and J. Shaham, _Effects of neutron-star superconductivity on magnetic monopoles and core field decay_ , _Phys. Rev. D_ 33 (1986) 2084.
  * (66) D.M. Sedrakian, A.D. Sedrakian and G.F. Zharkov, _Type i superconductivity of protons in neutron stars_ , _Monthly Notices of the Royal Astronomical Society_ 290 (1997) 203–207.
  * (67) B. Haskell and A. Sedrakian, _Superfluidity and superconductivity in neutron stars_ , _Astrophysics and Space Science Library_ (2018) 401–454.
  * (68) P. Landry, R. Essick and K. Chatziioannou, _Nonparametric constraints on neutron star matter with existing and upcoming gravitational wave and pulsar observations_ , _Physical Review D_ 101 (2020) .
  * (69) E. Berti. Private Communication, 2020.
  * (70) H.O. Silva, H. Sotani and E. Berti, _Low-mass neutron stars: universal relations, the nuclear symmetry energy and gravitational radiation_ , _Monthly Notices of the Royal Astronomical Society_ 459 (2016) 4378–4388.
  * (71) F. Özel and P. Freire, _Masses, radii, and the equation of state of neutron stars_ , _Annual Review of Astronomy and Astrophysics_ 54 (2016) 401 [https://doi.org/10.1146/annurev-astro-081915-023322].
  * (72) Super-Kamiokande collaboration, _Search for GUT monopoles at Super-Kamiokande_ , _Nucl. Phys. B Proc. Suppl._ 229-232 (2012) 540.
  * (73) K. Abe, Y. Hayato, T. Iida, M. Ikeda, C. Ishihara, K. Iyogi et al., _Solar neutrino results in super-kamiokande-iii_ , _Physical Review D_ 83 (2011) .
  * (74) G.W. Simon and R.B. Leighton, _Velocity Fields in the Solar Atmosphere. III. Large-Scale Motions, the Chromospheric Network, and Magnetic Fields._ , _The Astrophysical Journal_ 140 (1964) 1120.
  * (75) Z. Abraham and J. Iben, Icko, _More Solar Models and Neutrino Fluxes_ , _The Astrophysical Journal_ 170 (1971) 157.
  * (76) R. Arlt, A. Sule and G. Rüdiger, _Stability of toroidal magnetic fields in the solar tachocline_ , _Astronomy and Astrophysics_ 461 (2007) 295.
  * (77) MACRO collaboration, _Search for nucleon decays induced by GUT magnetic monopoles with the MACRO experiment_ , _Eur. Phys. J. C_ 26 163 [hep-ex/0207024].
  * (78) B.T. Draine, _Physics of the Interstellar and Intergalactic Medium_ (2011).
  * (79) F. "Timmes, “coldwd.tbz.”
  * (80) S. Grunau and V. Kagramanova, _Geodesics of electrically and magnetically charged test particles in the reissner-nordström space-time: Analytical solutions_ , _Physical Review D_ 83 (2011) .
  * (81) B. Carr, K. Kohri, Y. Sendouda and J. Yokoyama, _Constraints on primordial black holes_ , 2020.
  * (82) LIGO Scientific, Virgo collaboration, _GW190814: Gravitational Waves from the Coalescence of a 23 Solar Mass Black Hole with a 2.6 Solar Mass Compact Object_ , _Astrophys. J. Lett._ 896 (2020) L44 [2006.12611].
  * (83) B. Dasgupta, R. Laha and A. Ray, _Low mass black holes from dark core collapse_ , 2020.
  * (84) C.B. Richards, T.W. Baumgarte and S.L. Shapiro, _Relativistic bondi accretion for stiff equations of state_ , 2021.
  * (85) K. Jedamzik and L. Pogosian, _Relieving the hubble tension with primordial magnetic fields_ , 2020.
  * (86) R. Beck and R. Wielebinski, _Magnetic Fields in the Milky Way and in Galaxies_ , pp. 641–723 (2013), DOI [1302.5663].
  * (87) A.-C. Davis, M. Lilley and O. Törnkvist, _Relaxing the bounds on primordial magnetic seed fields_ , _Phys. Rev. D_ 60 (1999) 021301.
  * (88) L. Liu, O. Christiansen, Z.-K. Guo, R.-G. Cai and S.P. Kim, _Gravitational and Electromagnetic Radiations from Binary Black Holes with Electric and Magnetic Charges: I. Circular Orbits on a Poincaré Cone_ , 2008.02326.
  * (89) L. Liu, Z.-K. Guo, R.-G. Cai and S.P. Kim, _Merger rate distribution of primordial black hole binaries with electric charges_ , _Physical Review D_ 102 (2020) .
  * (90) L. Liu, Øyvind Christiansen, Z.-K. Guo, R.-G. Cai and S.P. Kim, _Gravitational and electromagnetic radiation from binary black holes with electric and magnetic charges: Elliptical orbits on a cone_ , 2020.
  * (91) P.C. Peters and J. Mathews, _Gravitational radiation from point masses in a keplerian orbit_ , _Phys. Rev._ 131 (1963) 435.
  * (92) M. Sasaki, T. Suyama, T. Tanaka and S. Yokoyama, _Primordial black holes—perspectives in gravitational wave astronomy_ , _Class. Quant. Grav._ 35 (2018) 063001 [1801.05235].
  * (93) T. Nakamura, M. Sasaki, T. Tanaka and K.S. Thorne, _Gravitational waves from coalescing black hole MACHO binaries_ , _Astrophys. J. Lett._ 487 (1997) L139 [astro-ph/9708060].