aainstitutetext: Tsung-Dao Lee Institute, Shanghai Jiao Tong University,
Shanghai 200240, Chinabbinstitutetext: Institute of Theoretical Physics,
Faculty of Physics, University of Warsaw, ul. Pasteura 5, 02-093 Warsaw,
Polandccinstitutetext: School of Physics & Astronomy, University of
Southampton, Southampton SO17 1BJ, UKddinstitutetext: Instituto de Física
Corpuscular, Universidad de Valencia and CSIC, Edificio Institutos
Investigación, C/Catedrático José Beltrán 2, 46980 Paterna, Spain

# Type-I two-Higgs-doublet model and gravitational waves from domain walls
bounded by strings

Bowen Fu b Anish Ghoshal c Stephen F. King d and Moinul Hossain Rahat
<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

The spontaneous breaking of a $U(1)$ symmetry via an intermediate discrete
symmetry may yield a hybrid topological defect of _domain walls bounded by
cosmic strings_. The decay of this defect network leads to a unique
gravitational wave signal spanning many orders in observable frequencies, that
can be distinguished from signals generated by other sources. We investigate
the production of gravitational waves from this mechanism in the context of
the type-I two-Higgs-doublet model extended by a $U(1)_{R}$ symmetry, that
simultaneously accommodates the seesaw mechanism, anomaly cancellation, and
eliminates flavour-changing neutral currents. The gravitational wave spectrum
produced by the string-bounded-wall network can be detected for $U(1)_{R}$
breaking scale from $10^{12}$ to $10^{15}$ GeV in forthcoming interferometers
including LISA and Einstein Telescope, with a distinctive $f^{3}$ slope and
inflexion in the frequency range between microhertz and hertz.

##  1 Introduction

After the discovery of the Higgs boson [1, 2], one of the next goals for the
collider experiments is to investigate the properties of the Higgs boson. The
2012 data has already shown an excess at 95 GeV [3] and next generation
collider experiments have been proposed to measure the Higgs boson
interactions precisely [4, 5, 6]. While revealing the nature of the standard
model (SM) Higgs boson, indication of new physics beyond the standard model
(BSM) may also appear.

One of the most attractive models that can be tested by measuring the
properties of the Higgs boson is the two-Higgs-doublet model (2HDM) [7]. In
2HDM, a serious potential problem is the tree-level flavour-changing neutral
currents (FCNCs). The necessary and sufficient condition to avoid FCNC is that
all fermions with the same quantum numbers couple to the same Higgs doublet
[8, 9], leading to limited types of the 2HDM [10]. A simple way to realise a
specific type of 2HDM is to impose an extra symmetry of the model. The most
popular possibility is to impose a $Z_{2}$ discrete symmetry, which can be
used to construct many different types of 2HDM [11]. However, as the SM is a
gauged theory, it is more interesting and technically elegant to consider a
gauge extension of the model symmetry.

In this work, we consider a gauged $U(1)$ extension of the SM gauge group as
the origin of natural flavour conservation [12, 13, 14, 15]. Furthermore, in
order to explain the neutrino oscillation, we include right-handed neutrinos,
which are responsible for generating the light neutrino masses through the
type-I seesaw mechanism [16, 17, 18, 19, 20]. Assuming the SM fermion content
is only extended by three right-handed neutrinos111Other types of 2HDM can be
implemented by adding vector-like fermions, see Ref. [12]., in the case of the
type-I 2HDM (but not other types), it is possible to control FCNCs by an
anomaly free gauged $U(1)$ symmetry. In general, type-I 2HDMs have a rich lab-
based phenomenology [21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34,
35, 36]. Here we shall focus on the implications of the gauged $U(1)$ type-I
2HDM for gravitational wave signatures.

In addition to the Higgs doublets, a complex scalar singlet is also added to
the scalar sector, which gives the right-handed neutrinos Majorana masses at
the seesaw scale. The vacuum expectation value (VEV) of the complex scalar
singlet breaks the $U(1)_{R}$ symmetry, under which only the right-handed
fermions are charged, to a residual $Z_{2}$ symmetry. Later, the $Z_{2}$
symmetry is further broken by the VEVs of Higgs doublets at the electroweak
(EW) scale. Motivated by experimental constraints [37, 38], we consider the
so-called “alignment limit”, where one of the scalar mass eigenstates recovers
the properties of the SM Higgs boson.

A natural outcome of this model is that the $U(1)_{R}$ symmetry is
spontaneously broken to a residual $Z_{2}$ symmetry at a high scale, motivated
by requiring heavy Majorana mass of the RHNs. This leads to the creation of a
cosmic string network. This network forms closed loops, and the vibrating
loops loose their energy via emitting gravitational waves. Later, when the
Higgs doublets get a VEV around the EW scale, the $Z_{2}$ symmetry is
spontaneously broken, leading to the creation of domain walls. These domain
walls are bounded by the pre-existing cosmic string network, creating a hybrid
network of “domain walls bounded by cosmic strings”. This topological defect
is unstable. Essentially, the cosmic string network loses energy by emitting
gravitational waves. The interaction between the surface tension of the wall
on the oscillating string loop accelerates its decay, resulting in an earlier
collapse of the string network. The net observational effect is the
modification of the gravitational wave spectrum at lower frequencies that
correspond to the late stage evolution of the defect network, while the high
frequency spectrum remain flat, as is also seen in the case of pure strings in
the absence of any intermediate symmetry breaking. In particular, the hybrid
network features a signal with $f^{3}$ slope in the infrared frequencies,
quite distinguishable from the pure string case. The amplitude of the high
frequency flat spectrum conveys information about the $U(1)_{R}$ breaking
scale, while the transition frequency where the infrared tail starts depends
on the intermediate $Z_{2}$ breaking scale. Since the latter is quite
constrained around the EW scale in the present model, we have a sharp
prediction for GW signals with distinct spectral slopes in the microHz – Hz
band.

The paper is structured as follows. In Sec. 2, we discuss the anomaly-free
$U(1)$ extension of the 2HDM with right-handed neutrinos. Section 3
investigates the emergence of the hybrid topological defect, and the
associated GW signal in this model. We report our main results in Sec. 4 and
finally conclude in Sec. 5.

##  2 Type-I 2HDM extended with a gauged $U(1)$ symmetry

We consider a gauged $U(1)$ extension of the SM gauge group in 2HDM with three
right-handed neutrinos. As a fundamental principle of extending the gauge
group, the chiral anomaly has to be cancelled. Suppose the charge of a
particle $i$ under the extra $U(1)$ symmetry is $X_{i}$. To achieve anomaly
cancellation, the charges of other fermions have to satisfy [12, 13, 14, 15]:

$\displaystyle
X_{q}=\frac{1}{2}(X_{d}+X_{u})\,,\,X_{l}=-\frac{3}{2}(X_{d}+X_{u})\,,\,X_{e}=-X_{d}-2X_{u}\,,\,X_{n}=-2X_{d}-X_{u}\,.$
(1)

$X_{q}$, $X_{u}$, $X_{d}$, $X_{l}$, $X_{e}$ and $X_{n}$ are the $U(1)$ charges
of left-handed quark doublet, right-handed up-type quark, right-handed down-
type quark, left-handed lepton doublet, right-handed charged lepton and right-
handed neutrino, respectively. Such a solution automatically leads to
$X_{u}-X_{q}=X_{q}-X_{d}=X_{l}-X_{e}=X_{n}-X_{l}$, meaning that the $U(1)$
charge of the Higgs doublet coupling to the fermions is enforced to be
$(X_{u}-X_{d})/2$. In order to avoid FCNC, the Higgs doublets must have
different quantum number, which ensures that only one of the Higgs doublets,
$\Phi_{2}$, can couple to the fermions while the other, $\Phi_{1}$, is
fermiophobic. The interaction between the fermiophobic Higgs doublet and the
fermions is forbidden as long as its charge is different from
$(X_{u}-X_{d})/2$. The type of 2HDM with such a structure in the Yukawa sector
is the well-known type-I 2HDM.

The $U(1)$ charge of the right-handed neutrinos forbids their Majorana masses.
In order to realise the seesaw mechanism, a scalar singlet $\phi$ is included
in the model, which couples to the right-handed neutrinos through
$\phi\overline{N^{c}}N$. The VEV of $\phi$ breaks the extra $U(1)$ symmetry
and gives mass to the right-handed neutrinos. As the charge of $\Phi_{1}$ is
not constrained due to the absence of its interaction with fermions, it is
possible for the scalar doublets to couple to the scalar singlet through
interaction of the type $(\Phi_{1}^{\dagger}\Phi_{2}\phi+\mathrm{h.c})$. Such
kind of possibility has been discussed in Ref. [15] with many concrete
examples.

In general, when the extra $U(1)$ symmetry is broken by the VEV of $\phi$,
there might be a residual discrete symmetry. The exact type of the residual
symmetry is determined by the charge of $\phi$ and the minimal charge unit of
the fields. For example, in the case that $X_{u}=1$ and $X_{d}=0$, the $U(1)$
charge of $\phi$ is 2 while the minimal charge unit of the fields is $1/2$.
After $\phi$ gets a VEV, the $U(1)$ symmetry would be broken into a $Z_{4}$.
Here, as the gravitational wave phenomena associated with the $Z_{2}$ domain
wall are the most well-understood, we choose $X_{u}=1$ and $X_{d}=-1$ for
simplicity. In such a case, the right-handed fermions are charged under the
$U(1)$ while the left-handed fermions are not, hence we will call this
$U(1)_{R}$.222Not to be confused with the continuous $R$-symmetry in the
context of supersymmetry, or the right-handed $U(1)$ symmetry arising from
left-right symmetric models. The particle content is listed in Table 1. After
$\phi$ gains a VEV, the $U(1)_{R}$ is broken into a $Z_{2}$ symmetry, which is
the simplest cyclic group

$\displaystyle
U(1)_{R}\stackrel{{\scriptstyle\langle\phi\rangle}}{{\longrightarrow}}Z_{2}.$
(2)

| ${u_{R}}_{\beta}$ | ${d_{R}}_{\beta}$ | ${Q}_{\alpha}$ | ${L}_{\alpha}$ | ${e_{R}}_{\beta}$ | $N_{R\beta}$ | $\Phi_{2}$ | $\Phi_{1}$ | $\phi$
---|---|---|---|---|---|---|---|---|---
$SU(2)_{L}$ | 1 | 1 | 2 | 2 | 1 | 1 | 2 | 2 | 1
$U(1)_{Y}$ | $\frac{2}{3}$ | $-\frac{1}{3}$ | $\frac{1}{6}$ | $-\frac{1}{2}$ | $-1$ | 0 | $\frac{1}{2}$ | $\frac{1}{2}$ | 0
$U(1)_{R}$ | $1$ | $-1$ | $0$ | $0$ | $-1$ | $1$ | $1$ | $-1$ | $-2$
residual $Z_{2}$ | $-$ | $-$ | $+$ | $+$ | $-$ | $-$ | $-$ | $-$ | $+$
Table 1:  Field content of the type-I 2HDM with gauged $U(1)_{R}$ symmetry and
its residual $Z_{2}$ subgroup.

With the particle content in Table 1, the allowed fermion Yukawa interactions
are

$\displaystyle\mathcal{L}_{Y}\supset
Y_{u}\overline{Q}\tilde{\Phi}_{2}u_{R}+Y_{d}\overline{Q}\Phi_{2}d_{R}+Y_{e}\overline{L}\Phi_{2}e_{R}+Y_{N}\overline{L}\tilde{\Phi}_{2}N_{R}+y_{N}\phi
N_{R}N_{R}+\text{h.c.},$ (3)

where $\tilde{\Phi}_{2}=-i\sigma_{2}\Phi_{2}^{*}$. Between the Higgs doublets,
only $\Phi_{2}$ is allowed to couple to the fermions. The most general scalar
potential allowed by the $SU(2)_{L}\times U(1)_{Y}\times U(1)_{R}$ symmetry is

$\displaystyle V(\Phi_{1},\,\Phi_{2},\,\phi)$ $\displaystyle=$ $\displaystyle
m_{11}^{2}\Phi_{1}^{\dagger}\Phi_{1}+m_{22}^{2}\Phi_{2}^{\dagger}\Phi_{2}+m_{\phi}^{2}|\phi|^{2}-m_{12}\left(\Phi_{1}^{\dagger}\Phi_{2}\phi+\mathrm{h.c}\right)$
(4)
$\displaystyle+\frac{\lambda_{1}}{2}(\Phi_{1}^{\dagger}\Phi_{1})^{2}+\frac{\lambda_{2}}{2}(\Phi_{2}^{\dagger}\Phi_{2})^{2}+\lambda_{3}(\Phi_{1}^{\dagger}\Phi_{1})(\Phi_{2}^{\dagger}\Phi_{2})+\lambda_{4}(\Phi_{1}^{\dagger}\Phi_{2})(\Phi_{2}^{\dagger}\Phi_{1})$
$\displaystyle+\lambda_{1\phi}(\Phi_{1}^{\dagger}\Phi_{1})|\phi|^{2}+\lambda_{2\phi}(\Phi_{2}^{\dagger}\Phi_{2})|\phi|^{2}+\lambda_{\phi}|\phi|^{4}\,.$

The potential is CP conserving as the phase of $m_{12}$ can be eliminated
through rephasing the scalar fields and all the other parameters are required
to be real by hermiticity. The $m_{12}$ term breaks the normal $B-L$ symmetry
explicitly and thus there is no spontaneous $U(1)_{B-L}$ symmetry
breaking.333Different from Ref. [15], the $B-L$ symmetry refers to the one
under which both of the Higgs doublets are uncharged. Despite the absence of a
strict $B-L$ symmetry, the Majorana mass of the right-handed neutrinos is
still protected by the $U(1)_{R}$ symmetry. After $\phi$ gets a VEV
$\langle\phi\rangle=v_{\rm{M}}$, the potential of the Higgs doublets becomes

$\displaystyle V(\Phi_{1},\,\Phi_{2})$ $\displaystyle=$
$\displaystyle\tilde{m}_{11}^{2}\Phi_{1}^{\dagger}\Phi_{1}+\tilde{m}_{22}^{2}\Phi_{2}^{\dagger}\Phi_{2}-\tilde{m}_{12}^{2}\left(\Phi_{1}^{\dagger}\Phi_{2}+\mathrm{h.c}\right)$
(5)
$\displaystyle+\frac{\lambda_{1}}{2}(\Phi_{1}^{\dagger}\Phi_{1})^{2}+\frac{\lambda_{2}}{2}(\Phi_{2}^{\dagger}\Phi_{2})^{2}+\lambda_{3}(\Phi_{1}^{\dagger}\Phi_{1})(\Phi_{2}^{\dagger}\Phi_{2})+\lambda_{4}(\Phi_{1}^{\dagger}\Phi_{2})(\Phi_{2}^{\dagger}\Phi_{1})\,,$

where $\tilde{m}_{11}^{2}=m_{11}^{2}+\lambda_{1\phi}v_{\rm{M}}^{2}$,
$\tilde{m}_{22}^{2}=m_{22}^{2}+\lambda_{2\phi}v_{\rm{M}}^{2}$ and
$\tilde{m}_{12}^{2}=m_{12}v_{\rm{M}}$. Unlike the most general CP-conserving
potential for 2HDMs, the interaction
$\frac{\lambda_{5}}{2}((\Phi_{1}^{\dagger}\Phi_{2})^{2}+\mathrm{h.c})$ is
still absent due to the $U(1)_{R}$ symmetry.

The most general 2HDM admits CP-breaking vacua, where the VEVs of the Higgs
doublets have a relative phase. However, as all the parameters in the scalar
potential in Eq. (5) are real, the vacua are always real and can be expressed
as $\langle\Phi_{i}\rangle=v_{i}/\sqrt{2}$, where
$\sqrt{v_{1}^{2}+v_{2}^{2}}=v_{\rm{SM}}\sim 246$ GeV.444A relative phase
$\varphi$ only affects the potential by a factor $\cos\varphi$ in the
$\tilde{m}_{12}^{2}$ term. Since $\tilde{m}_{12}^{2}$ is real, minimising the
potential leads to $\cos\varphi=1$ for $\tilde{m}_{12}^{2}>0$ and
$\cos\varphi=-1$ for $\tilde{m}_{12}^{2}<0$.  The ratio $v_{2}/v_{1}$ is
defined as $\tan\beta$. Rotating the Higgs fields by an angle $\beta$, we
obtain the so-called “Higgs basis”, where only one of the doublets gets a
nonzero VEV. Another popular basis for 2HDM is the “mass eigenbasis”, found by
a rotating the fields by an angle $\alpha$. However, the properties of SM
Higgs are not recovered by any field in either of these choices. Instead, a
combination defined as

$\displaystyle H_{\rm SM}=H\cos{(\beta-\alpha)}+h\sin{(\beta-\alpha)}$ (6)

recovers both the gauge interaction and Yukawa interaction of the SM Higgs
boson, where $H=(\Phi_{1}\cos\alpha+\Phi_{2}\sin\alpha)$ and
$h=(\Phi_{2}\cos\alpha-\Phi_{1}\sin\alpha)$ are the two mass eigenstates. For
details see Appendix A.

###  2.1 The alignment limit and the domain wall solution

As the properties of the Higgs boson discovered in 2012 are consistent with
the predication of the SM at around $2\sigma$ level [37, 38], we first
consider the “alignment limit” [39, 40], where one of the Higgs mass
eigenstates recovers the properties of the SM Higgs boson. This is achieved,
for example, by imposing the relation $\alpha=\beta-\pi/2$, which makes $h$
the SM-like Higgs boson, i.e. $h=H_{\rm SM}$ in Eq. (6).

The scalar potential is symmetric under the surviving $Z_{2}$ in Eq. (2) for
$\Phi_{i}\to-\Phi_{i}$. When the Higgs doublets get VEV $v_{i}$, this symmetry
is spontaneously broken, resulting in regions of space dubbed as “domains”
randomly adopting the $+v_{i}$ or $-v_{i}$ vacua. These regions are separated
by walls with a static planar configuration, taken to be perpendicular to the
$z$ direction without loss of generality. The field equation for the domain
walls can be expressed as [41]

$\displaystyle\frac{d^{2}\Phi_{i}}{dz^{2}}=\frac{\partial
V}{\partial\Phi_{i}}\,,\quad i=1,2\,,$ (7)

with the boundary conditions $\Phi_{i}(z\to\pm\infty)=\pm v_{i}$ at two
different vacua. For a general potential with multiple scalar fields, one
needs to consider the parallel and perpendicular parts of the equations of
motion, which are

$\displaystyle\frac{d^{2}x}{dz^{2}}$ $\displaystyle=$
$\displaystyle\frac{dV[\vec{\Phi}(x)]}{dx}\,,$ (8)
$\displaystyle\frac{d^{2}\vec{\Phi}}{dx^{2}}\left(\frac{dx}{dz}\right)^{2}$
$\displaystyle=$ $\displaystyle\nabla_{\perp}V(\vec{\Phi})\,,$ (9)

and apply the path deformation method [42, 35] to find the solution. However,
in the alignment limit, the correct path is simply the straight line
connecting the two minima in the field space (the blue path in Fig. 1a).
Indeed, once we choose the path to be $\Phi_{1}=\cos\beta x$ and
$\Phi_{2}=\sin\beta x$, Eq. (9) is automatically satisfied regardless of the
relation between $x$ and $z$ (see Appendix B for details). Then the problem is
reduced to finding the domain wall solution for Eq. (8), which is simply

$\displaystyle\frac{d^{2}x}{dz^{2}}$ $\displaystyle=$
$\displaystyle\frac{m_{h}^{2}}{2v_{\rm SM}^{2}}x(x^{2}-v_{\rm SM}^{2})=\lambda
x(x^{2}-v_{\rm SM}^{2})\,$ (10)

with $m_{h}=125$ GeV the mass of the SM-like Higgs $h$ and $\lambda\equiv
m_{h}^{2}/2v_{\rm SM}^{2}$ the SM Higgs quartic coupling. Here we have taken
$v_{2}\sim v_{\rm{SM}}$. Eq. (10) admits a well-known analytical solution [43]

$\displaystyle x(z)=v_{\rm SM}\tanh{\left(z/\Delta\right)}\,,$ (11)

where $\Delta=\sqrt{2/\lambda}/v_{\rm SM}$ is the thickness of the wall. As
usual, the surface energy density of the wall is

$\displaystyle\mathcal{E}=\frac{2}{3}\sqrt{2\lambda}v_{\rm SM}^{3}\,,$ (12)

which is around $5.04\times 10^{6}\,\text{GeV}^{3}$ given the SM Higgs VEV and
self coupling.

(a)

(b)

Figure 1: Path in the field space (left) and energy density across the domain
wall (right). The energy density $\mathcal{E}$ for $\delta=+0.2$ and
$\delta=-0.1$ are $5.75\times 10^{6}\,\text{GeV}^{3}$ and $5.25\times
10^{6}\,\text{GeV}^{3}$, respectively.

###  2.2 Deviation from the alignment limit

Although the experiments indicate the alignment limit is a good approximation,
a deviation has not been ruled out yet. Allowing a deviation from the
alignment limit, the domain wall solution would no longer follow the straight
line between the two equivalent minima (see Appendix B for details). The
domain wall solution can be found by the path deformation method (see for
example [42, 35]), or alternatively, by the overshoot/undershoot method (see
for example [44, 42]). As the scalar potential possesses a central symmetry,
the domain wall path has to pass the center, which is the origin. The
magnitude of the field profile derivative can be determined everywhere in the
field space as $|d\vec{\phi}/dz|=\sqrt{2\Delta V}$ in the spirit of the
overshoot/undershoot method. Eventually the problem is reduced to finding the
direction of the path at the origin. The path exactly connect the two minima
for the correct direction, while the wrong direction never pass the minima
without any ends, as shown in Fig. 2.

Figure 2: Paths found by the overshoot/undershoot method. The blue path
exactly connects the two minima (green points). The red and purple lines show
two paths with the wrong direction. The yellow point is the symmetry center.

Using both methods, we solve the coupled system of field equations (8) and
(9), and explore how much a deviation from the SM-like Higgs limit can change
the domain wall energy density.

An illustrative example of this is shown in Fig. 1a. After the spontaneous
symmetry breaking, there are five physical scalar fields, namely, two neutral
scalars (one of which, $h$, aligns with the SM Higgs in the alignment limit),
a charged scalar, and a pseudoscalar. For the present example, we choose a
benchmark point where the masses of the neutral scalar different from the SM-
like Higgs, the charged scalar and the pseudoscalar are $410$ GeV, $380$ GeV
and $400$ GeV, respectively, while $\tan\beta=5$. Details on the relations
between the masses of physical scalar fields and model parameters are
discussed in Appendix A. The path in the field space for the domain wall
solution follows a straight line in the perfect alignment limit, as shown in
Fig. 1a. On the other hand, allowing a deviation $\delta$ from the alignment
limit, i.e., $\alpha=\beta-\pi/2+\delta$, the path changes. In particular, we
consider two extremes of the deviation angle $\delta=+0.2$ and $\delta=-0.1$
which would reduce the coupling between $h$ and gauge bosons by around 2% and
0.5%, respectively.555A deviation $\delta$ from the alignment limit changes
both the gauge interaction ($hW^{\mu}W_{\mu}$) and the Yukawa interaction
($h\overline{t}t$) of $h$ by a factor $\cos\delta$. $\delta=+0.2$ and
$\delta=-0.1$ correspond to $\cos\delta\simeq 0.98$ and $\cos\delta\simeq
0.995$, respectively. Such deviations represent the limit allowed by the the
experimental constraints [38, 37] while preserving the stability of the
vacua.666A deviation $\delta=-0.2$ leads to metastable vacua for the benchmark
case. The energy density in the direction perpendicular to the wall is
presented in Fig. 1b. In the presence of the deviation, the energy at the
center of the domain wall becomes larger, leading to an increase in the total
surface energy density. However, we find that the variation of the domain wall
surface energy density is only about $15\%$, leading to a negligible effect in
the gravitational wave phenomenology.

##  3 Domain walls bounded by cosmic strings

In the type-I 2HDM extended by a gauged $U(1)_{R}$, the symmetry breaking
pattern is as follows. The scalar $\phi$ gets a nonzero VEV at a high scale
and gives large mass to the right-handed neutrino. This spontaneously breaks
the $U(1)_{R}$ symmetry and creates a cosmic string network [45]. We assume
that this happens in a radiation dominated Universe after inflation has ended.
This ensures that the string network is not diluted by inflation. 777Unlike
other topological defects, string networks are not completely diluted away by
inflation, and can even regrow to generate observable gravitational wave
signatures. See [46, 47, 48] for details. In case that the string is
completely diluted, hybrid topological defects can still appear through
nucleation of holes bounded by strings on the domain wall [49, 50]. However,
such topological defects can only lead to observable gravitational wave for
quasi-degenerate string scale and wall scale, with inflation scale in the
middle.

Because of the $U(1)_{R}$ charge assignment of the fields, as given in Table
1, the model has a residual $Z_{2}$ symmetry at the $\phi$ vacuum, under which
both Higgs doublets $\Phi_{1}$ and $\Phi_{2}$ are odd. When either of them
gets a nonzero VEV, the $Z_{2}$ symmetry is spontaneously broken, leading to
the creation of a domain wall network [45].

While the creation of cosmic strings and domain walls from their respective
symmetry breakings have been widely studied in the context of various BSM
models, the situation above is somewhat different. Here the $Z_{2}$ is not a
symmetry of the original Lagrangian, but is a residual symmetry that occurs at
an intermediate stage in the symmetry breaking chain. When the $U(1)_{R}$
symmetry is broken, the VEV of the scalar field with $U(1)_{R}$ charge 2 can
have a phase $\theta(x)$ that can vary between $0$ and $\pi$ with the
spacetime. The same applies to the other fields in the theory. After the
$U(1)_{R}$ symmetry breaking, the phase changes from $\theta(x)$ to
$\theta(x)+\pi$ around a string, as shown in Fig. LABEL:fig2a. However, those
fields with $U(1)_{R}$ charge 1 are not invariant when $\theta(x)$ changes to
$\theta(x)+\pi$. Instead, there is a sign change in those fields, leading to a
residual $Z_{2}$ symmetry. As such, the created domain walls after $Z_{2}$
breaking are not infinite planar walls, rather circular walls that fill in the
space between string loops, as shown in Fig. LABEL:fig2b.888We focus on string
loops. In case that the wall is formed between long strings, the wall would
promote the loop formation [51]. This creates a hybrid topological defect
dubbed as “walls bounded by strings” [52].

(a)

(b)

Figure 3: Illustration of the appearance of the hybrid topological defect. The
spontaneous breaking of $U(1)_{R}$ into a residual $Z_{2}$ symmetry creates a
string network. When the $Z_{2}$ symmetry is spontaneously broken later,
domain walls appear which fill up the space between the string loops.

A detailed description of the evolution of this hybrid topological defect can
be found in Ref. [53]. Here we will briefly discuss the relevant physics and
the observational signatures of the defect.

A useful parameter that dictates the evolution of the defect is defined as
$R_{c}\equiv\mu/\mathcal{E}$, where $\mu\sim v_{\rm M}^{2}$ is the string
tension, $v_{\rm M}$ being the $U(1)_{R}$ breaking scale, and $\mathcal{E}\sim
v_{\rm{SM}}^{3}$ is the surface energy density of the domain wall,
$v_{\rm{SM}}$ being the $Z_{2}$ breaking scale, and assuming a
$\mathcal{O}(1)$ coupling constant of the scalar field breaking $Z_{2}$. After
the $U(1)_{R}$ symmetry breaking, horizon size long strings are created that
intersect and form string loops. These circular loops oscillate non-
relativistically and lose their energy by emitting gravitational waves. When
domain walls are formed at $t_{\rm DW}\approx M_{\rm Pl}C/v_{\rm{SM}}^{2}$,
where $C=(8\pi^{3}g_{\star}/90)^{-1/2}$, and fill up the space between the
string loops, they tend to make the string motion ultra-relativistic as long
as the loop size is above $R_{c}$. This makes the wall-string network collapse
earlier than the case of pure string loops. If $t_{\rm DW}<R_{c}$, walls do
not dominate the string dynamics immediately after creation. As such, the
hybrid network becomes ultra-relativistic and collapses after some time,
$t_{\star}\sim R_{c}$. On the other hand, if $t_{\rm DW}>R_{c}$, the hybrid
network becomes ultra-relativistic as soon as the walls are created and the
network collapses at $t_{\star}\sim t_{\rm DW}$. We can, therefore, define the
time for the collapse of the hybrid network as
$t_{\star}\equiv\max(R_{c},t_{\rm DW})$ [54].

For a loop of radius $R$ with circular length $l=2\pi R$, the rate of energy
loss is given by

$\displaystyle\frac{dE}{dt}=-\Gamma(l)G\mu^{2},$ (13)

where $G$ is the Newton’s constant. The function
$\Gamma(l)\approx\Gamma_{s}=50$ when $l\ll 2\pi R_{c}$ in the pure string
limit, and $\Gamma(l)=3.7(l/(2\pi R_{c}))^{2}$ for $l\gg 2\pi R_{c}$ in the
pure wall limit. The behavior around $l\sim 2\pi R_{c}$ can be approximated
using a smooth interpolation function, as shown in Ref. [53].

Because of the energy loss, loops forming at time $t_{k}$ with an initial size
$l_{k}=\alpha t_{k}$ slowly decrease in size with time. Here $\alpha\approx
0.1$ is the ratio between loop formation length and horizon size, and its
value is determined from simulations [55, 56]. For the domain walls bounded by
strings, the size of a string loop emitting gravitational waves at time
$\tilde{t}$ is $\tilde{l}=2\pi\tilde{R}$, and we have

$\displaystyle G\mu(\tilde{t}-t_{k})=\int_{\tilde{l}}^{\alpha
t_{k}}dl^{\prime}\frac{1+\frac{l^{\prime}}{2\pi R_{c}}}{\Gamma(l^{\prime})}.$
(14)

Here the integration variable $l^{\prime}$ can be interpreted as the
instantaneous length of the loop emitting gravitational wave at time
$t^{\prime}$, where $t_{k}<t^{\prime}<\tilde{t}$. We further assume that at
time $\tilde{t}$, the newly created loops are larger than the pre-existing
loops emitting gravitational wave at time $\tilde{t}$

$\displaystyle\alpha\tilde{t}>\tilde{l}.$ (15)

We can decompose the gravitational wave into Fourier modes $\tilde{f}_{k}=\xi
k/\tilde{l}$, where the parameter $\xi$ varies between $2$ in the pure string
limit ($\tilde{l}\ll 2\pi R_{c}$) and $\pi$ in the pure wall limit
$\tilde{l}\gg 2\pi R_{c}$ [53]. This mode is redshifted to the frequency
$f_{k}$ in today’s interferometers

$\displaystyle f_{k}=\frac{a(\tilde{t})}{a(t_{0})}\tilde{f}_{k},$ (16)

where $a(t)$ is the scale factor which is defined as $a(t_{0})\equiv 1$ where
$t_{0}$ is the present time.

The gravitational wave spectrum observed today is given by summing over all
Fourier modes

$\displaystyle\Omega_{\rm GW}h^{2}(f)$
$\displaystyle=\sum_{k}\frac{G\mu^{2}}{\rho_{c}}\int_{t_{\rm
sc}}^{t_{0}}d\tilde{t}\left(\frac{a(\tilde{t})}{a(t_{0})}\right)^{5}\left[\frac{\mathcal{F}C_{\rm
eff}(t_{k})}{\alpha
t_{k}^{4}}\left(\frac{a(t_{k})}{a(\tilde{t})}\right)^{3}\right]\mathcal{P}_{k}\frac{\xi
k}{f}\left[1+\frac{1}{2\pi R_{c}}\frac{\xi
k}{f}\frac{a(\tilde{t})}{a(t_{0})}\right]$
$\displaystyle\times\frac{\Gamma(\alpha
t_{k})\theta{(t_{\star}-t_{k})}}{\Gamma(\alpha
t_{k})G\mu+\alpha\left(1+\frac{\alpha t_{k}}{2\pi R_{c}}\right)}.$ (17)

Here $\rho_{c}=3H_{0}^{2}/(8\pi G)$ is the critical energy density today.
$t_{\rm sc}\sim 10^{-22}$ s is the time when the network reaches a scaling
regime where the energy density of the hybrid network becomes constant in the
expanding Universe. $\mathcal{F}\approx 0.1$ is the fraction of energy that is
transferred by the string network into loops [56]. $C_{\rm eff}$ is the loop
formation efficiency which equals $5.7$ during the radiation-dominated era,
and $0.5$ during the matter-dominated era. $\mathcal{P}_{k}=k^{-q}/\zeta(q)$
is the fractional power radiated by the $k$th mode of an oscillating string
loop with cusp, and $\zeta(x)=\sum_{n=1}^{\infty}n^{-x}$ is the Riemann zeta
function. The power spectral index is found to be $q=4/3$ for string loops
containing cusps [57, 58].

In Eq. (17), we are integrating over the gravitational wave modes emitted by
the hybrid network from an early time when the network goes to a scaling
regime, $t_{\rm sc}$, to today $t_{0}$. Here the integration variable
$\tilde{t}$ represents the time when a loop is radiating gravitational wave,
and $t_{k}$ is the time when this particular loop was created, which can be
obtained by solving Eq. (14). The Heaviside $\theta$ function
$\theta(t_{\star}-t_{k})$ ensures that we only sum over the contribution from
loops which were created before the network collapses at $t_{\star}$.

We have developed a public code package CosmicStringGW, available on GitHub at
/moinulrahat/CosmicStringGW, for calculating the gravitational wave spectrum
from pure and hybrid defect networks of cosmic strings.

##  4 Gravitational wave signature

We expect the $U(1)_{R}$ symmetry to be spontaneously broken at a high scale
above $10^{10}$ GeV, inspired by a high scale of right-handed neutrino masses,
and the $Z_{2}$ symmetry to be broken at the electroweak scale. Hence we focus
on keeping the $Z_{2}$ breaking scale fixed at $v_{\rm SM}=246$ GeV, while
varying the $U(1)_{R}$ breaking scale.999We will adopt the perfect alignment
limit when calculating the domain wall energy density, as we have verified
that the allowed deviation has unnoticeable impact on the gravitational wave
signals. The resulting gravitational waves from the string-wall network,
calculated from evaluating Eq. (17) up to $k=10^{5}$ terms, is shown in Fig. 4
with solid lines. In the same plot we also show the spectrum corresponding to
the pure string case, if the $U(1)_{R}$ symmetry were completely broken at the
high scale.101010This is achieved by taking the $R_{c}\to\infty$ limit in Eq.
(17). Projected sensitivity (SKA [59], $\mu$-Ares [60], LISA [61], DECIGO [62,
63], BBO [64], AEDGE [65], AION [66], CE [67], ET [68] and future upgrades of
LVK), upper bounds (LVK [69, 70, 71]) from various interferometers and recent
results from pulsar timing arrays NANOGrav [72, 73] and EPTA [74, 75] are
shown for comparison. Astrophysical foregrounds for the relevant frequency
range are not shown, and can be found in Ref. [76]. Signal-to-noise ratio
(SNR) at some of the detectors is tabulated in Appendix C.

Figure 4: GW spectrum from domain walls bounded by cosmic strings (solid
lines). Dotted lines represent the corresponding pure string induced GW
spectrum in the absence of an intermediate $Z_{2}$ breaking. See text for
details.

We notice the following key characteristics of the gravitational wave spectrum
generated from the hybrid network of domain walls bounded by cosmic strings.

  * •

Compared to the pure string case (dotted line), the signals from the hybrid
network (solid line) have an infrared (IR) tail at larger frequencies. This is
caused by the collapse of the network after creation of the domain walls that
drive the boundary strings to ultra-relativistic dynamics. The frequency where
the IR tail starts is determined by the collapse time $t_{\rm DW}$, which in
our model depends on $R_{c}\sim v_{\rm{M}}^{3}/v_{\rm SM}$ (since
$R_{c}>t_{\rm DW}$ for $v_{\rm{M}}\gtrsim 10^{10}$ GeV). Higher $R_{c}$
implies that the collapse happens later, allowing the flat part of the signal
to continue to lower frequencies (see Fig.4).

  * •

The spectral slope of the hybrid-defect signal is steeper compared to the pure
string case. This has important implications for interferometers. The steep
$f^{3}$ slope of the IR tail can be distinguished from the pure string case in
planned interferometer sensitivities in the $\mu$Hz to Hz frequencies.

  * •

The IR tail appearing at relatively higher frequencies compared to the pure
string case implies that non-observation at pulsar timing arrays (PTA) does
not rule out signals from the hybrid defect for $v_{\rm M}\gtrsim 10^{14}$
GeV. The upper bound $v_{\rm M}\lesssim 10^{15.5}$ GeV from LVK applies to
gravitational waves from both pure and hybrid cosmic string defects, as both
signals are identical at LVK frequencies.

  * •

While the pure string signal probes a single symmetry breaking scale [77, 78,
79, 80, 81], signals from hybrid defects like domain walls bounded by cosmic
strings can simultaneously probe two symmetry breaking scales ($U(1)_{R}$ and
$Z_{2}$). This offers a unique opportunity to test BSM scenarios where such
scales are separated by orders of magnitude.

  * •

Given that the $Z_{2}$ breaking scale in the current model is fixed at the
electroweak scale, the model predicts a signal with spectral slope $f^{3}$ in
the $\mu$Hz to Hz range for $U(1)_{R}$ breaking scale $10^{12}-10^{15}$ GeV.
Together with a flat spectrum at frequencies higher than the IR tail, this is
a unique signature of this model.

In summary, the hybrid network generated from the sequential breaking of the
$U(1)_{R}$ symmetry via the intermediate $Z_{2}$ symmetry results in a unique
GW signal distinguishable from other new physics signals in the observable
frequency range, and can be probed in upcoming GW interferometers.

##  5 Conclusion

We have investigated the emergence of a hybrid topological defect, namely
“domain walls bounded by cosmic strings”, and the unique gravitational wave
signal that emerges from it in the context of the type-I two-Higgs-doublet
model (2HDM) extended with a gauged $U(1)_{R}$ symmetry. The phenomenological
motivation for the $U(1)_{R}$ extension stems from the necessity to
incorporate neutrino masses and mixings via the seesaw mechanism. The model
simultaneously cancels the chiral anomaly and eliminates the flavour-changing
neutral currents. An interesting outcome of the $U(1)_{R}$ charge assignment,
motivated by spontaneous generation of the right-handed neutrino Majorana
masses, is that $U(1)_{R}$ is broken to an intermediate $Z_{2}$ symmetry,
which is further broken around the electroweak scale by the Higgs doublets.
Such a chain of symmetry breaking triggers a hybrid topological defect of
“domain walls bounded by cosmic strings”, created by the domain walls of
$Z_{2}$ breaking being encircled by pre-existing cosmic string loops
originated from the earlier $U(1)_{R}$ breaking.

We have studied the cosmological footprint of this model in terms of the decay
of the hybrid network of walls bounded by strings into gravitational waves.
The defect network evolves by losing energy via gravitation wave emission,
where the observable frequency of the signal is inversely proportional to the
time of evolution. The decay of the cosmic string loops created after
$U(1)_{R}$ breaking yields a flat spectrum of gravitational waves at high
frequencies. The additional tension exerted by the bounded domain walls
arising after $Z_{2}$ breaking causes an accelerated decay of the loops,
resulting in a characteristic $f^{3}$ slope of the signal in lower
frequencies. The pivot frequency where this spectral change occurs is
determined by the ratio of the two symmetry breaking scales, while the high
frequency flat amplitude is representative of the higher scale. In the present
model, the $Z_{2}$ is broken around the electroweak scale. This leads to a
prediction for the $f^{3}$ behaviour to be observed in the microHz to Hz band
with a flat spectrum in the ultraviolet, when the $U(1)_{R}$ breaking scale is
in between $10^{12}$ GeV to $10^{15}$ GeV. In this sense, our study provides a
way to test the model in upcoming gravitational wave interferometers such as
LISA and ET, complementing the collider results.

## Acknowledgement

We thank Yongcheng Wu for helpful discussions. MHR acknowledges financial
support from the Spanish AEI-MICINN
PID2020-113334GB-I00/AEI/10.13039/501100011033 and Generalitat Valenciana
project CIPROM/2021/054. SFK thanks CERN for hospitality and acknowledges the
STFC Consolidated Grant ST/L000296/1 and the European Union’s Horizon 2020
Research and Innovation programme under Marie Skłodowska-Curie grant agreement
HIDDeN European ITN project (H2020-MSCA-ITN-2019//860881-HIDDeN).

##  Appendix A The scalar sector of 2HDM

The minimisation conditions for the Higgs potential are given by

$\displaystyle\frac{\partial V}{\partial v_{1}}$
$\displaystyle=\tilde{m}_{11}^{2}v_{1}-\tilde{m}_{12}^{2}v_{2}+\frac{\lambda_{1}}{2}v_{1}^{3}+\frac{\lambda_{34}}{2}v_{1}v_{2}^{2}=0\,,$
(18a) $\displaystyle\frac{\partial V}{\partial v_{2}}$
$\displaystyle=\tilde{m}_{22}^{2}v_{2}-\tilde{m}_{12}^{2}v_{1}+\frac{\lambda_{2}}{2}v_{2}^{3}+\frac{\lambda_{34}}{2}v_{1}^{2}v_{2}=0\,,$
(18b)

where $\lambda_{34}=\lambda_{3}+\lambda_{4}$. After the spontaneous symmetry
breaking, there are five physical scalar fields: two neutral scalars ($H$ and
$h$), a charged scalar ($H^{\pm}$), and one pseudoscalar ($A$). The mass
spectrum of the scalar particles are

$\displaystyle m_{\pm}^{2}$
$\displaystyle=\tilde{m}_{12}^{2}(\tan\beta+\cot\beta)-\lambda_{4}v_{\rm
SM}^{2}\,,$ (19) $\displaystyle m_{A}^{2}$
$\displaystyle=\tilde{m}_{12}^{2}(\tan\beta+\cot\beta)\,,$ (20)
$\displaystyle\begin{pmatrix}m_{H}^{2}&0\\\ 0&m_{h}^{2}\end{pmatrix}$
$\displaystyle=\begin{pmatrix}c_{\alpha}&s_{\alpha}\\\
-s_{\alpha}&c_{\alpha}\end{pmatrix}\begin{pmatrix}\tilde{m}_{12}^{2}\tan\beta+\lambda_{1}v_{1}^{2}&-\tilde{m}_{12}^{2}+\lambda_{34}v_{1}v_{2}\\\
-\tilde{m}_{12}^{2}+\lambda_{34}v_{1}v_{2}&\tilde{m}_{12}^{2}\cot\beta+\lambda_{2}v_{2}^{2}\end{pmatrix}\begin{pmatrix}c_{\alpha}&-s_{\alpha}\\\
s_{\alpha}&c_{\alpha}\end{pmatrix}\,,$ (21)

where $c_{\alpha}=\cos\alpha$ and $s_{\alpha}=\sin\alpha$.

In total, there are 7 model parameters and 7 physical parameters in the Higgs
sector, as listed in Table 2.

Model parameters | Physical parameters
---|---
$\tilde{m}_{11,22,12},\,\lambda_{1,2,3,4}$ | $m_{\pm},\,m_{A},\,m_{H},\,m_{h},\,v_{\rm SM}\,,\beta,\,\alpha$
Table 2:  Parameters in the Higgs sector

The physical parameters and the model parameters are related by

$\displaystyle\tilde{m}_{11}^{2}$
$\displaystyle=\frac{1}{2}\sec\beta(-m_{H}^{2}\cos\alpha\cos(\alpha-\beta)-m_{h}^{2}\sin\alpha\sin(\alpha-\beta)+m_{A}^{2}\sin\beta\sin
2\beta)\,,$ (22) $\displaystyle\tilde{m}_{22}^{2}$
$\displaystyle=\frac{1}{2}\csc\beta(-m_{H}^{2}\sin\alpha\cos(\alpha-\beta)+m_{h}^{2}\cos\alpha\sin(\alpha-\beta)+m_{A}^{2}\cos\beta\sin
2\beta)\,,$ (23) $\displaystyle\tilde{m}_{12}^{2}$
$\displaystyle=\frac{m_{A}^{2}}{\tan\beta+\cot\beta}\,,$ (24)
$\displaystyle\lambda_{1}$ $\displaystyle=\frac{1}{v_{\rm
SM}^{2}\cos^{2}\beta}\left(m_{H}^{2}\cos^{2}\alpha+m_{h}^{2}\sin^{2}\alpha-
m_{A}^{2}\sin^{2}\beta\right)\,,$ (25) $\displaystyle\lambda_{2}$
$\displaystyle=\frac{1}{v_{\rm
SM}^{2}\sin^{2}\beta}\left(m_{H}^{2}\sin^{2}\alpha+m_{h}^{2}\cos^{2}\alpha-
m_{A}^{2}\cos^{2}\beta\right)\,,$ (26) $\displaystyle\lambda_{3}$
$\displaystyle=\frac{1}{v_{\rm SM}^{2}\sin 2\beta}(m_{H}^{2}-m_{h}^{2})\sin
2\alpha+\frac{m_{\pm}^{2}}{v_{\rm SM}^{2}}\,,$ (27) $\displaystyle\lambda_{4}$
$\displaystyle=\frac{m_{A}^{2}-m_{\pm}^{2}}{v_{\rm SM}^{2}}\,.$ (28)

In the alignment limit, these relations between model and physical parameters
are reduced to

$\displaystyle\tilde{m}_{11}^{2}=\frac{1}{2}(-m_{h}^{2}+2m_{A}^{2}\sin^{2}\beta)\,,$
$\displaystyle\tilde{m}_{22}^{2}=\frac{1}{2}(-m_{h}^{2}+2m_{A}^{2}\cos^{2}\beta)\,,$
$\displaystyle\quad\tilde{m}_{12}^{2}=\frac{m_{A}^{2}}{\tan\beta+\cot\beta}\,,$
$\displaystyle\lambda_{1}=\frac{(m_{H}^{2}-m_{A}^{2})\tan^{2}\beta+m_{h}^{2}}{v_{\rm
SM}^{2}}\,,$
$\displaystyle\lambda_{2}=\frac{(m_{H}^{2}-m_{A}^{2})\cot^{2}\beta+m_{h}^{2}}{v_{\rm
SM}^{2}}\,,$
$\displaystyle\lambda_{3}=\frac{m_{\pm}^{2}-m_{H}^{2}+m_{h}^{2}}{v_{\rm
SM}^{2}}\,,$ $\displaystyle\lambda_{4}=\frac{m_{A}^{2}-m_{\pm}^{2}}{v_{\rm
SM}^{2}}\,.$ (29)

The stability of the vacuum can be ensured by $\lambda_{1,2}>0$,
$\lambda_{3}>-\sqrt{\lambda_{1}\lambda_{2}}$ and
$\lambda_{3}+\lambda_{4}>-\sqrt{\lambda_{1}\lambda_{2}}$.

##  Appendix B Equations of motion for domain wall in the Higgs basis

In the 2HDM, there is a basis of the Higgs doublets under which only one of
the Higgs doublets obtains a VEV, which is known as the Higgs basis. The Higgs
basis can be obtained by the following rotation of the fields:

$\displaystyle\Phi_{1}\rightarrow S_{1}\cos\beta-
S_{2}\sin\beta\,,\quad\Phi_{2}\rightarrow S_{1}\sin\beta+S_{2}\cos\beta\,.$
(30)

In the Higgs basis, the scalar potential Eq. (5) would be minimised at
$(s_{1},s_{2})=(v_{\rm SM},0)$, where $s_{i}\equiv\sqrt{2}\langle
S_{i}\rangle$.

Under the alignment limit, the path of the field value for the domain wall
solution happens to be a straight line connecting the two $Z_{2}$ symmetric
vacua. Along the $S_{1}$ direction ($s_{2}=0$), the gradient of the potential
reads

$\displaystyle\left.\frac{\partial V}{\partial
s_{1}}\right|_{s_{2}=0}=\frac{m_{h}^{2}}{2v_{\rm
SM}^{2}}s_{1}(s_{1}^{2}-v_{\rm SM}^{2})\,,\quad\left.\frac{\partial
V}{\partial s_{2}}\right|_{s_{2}=0}=0\,,$ (31)

where $m_{h}$ is the mass of $h$. Therefore the equation of motion for
$s_{2}$, which is $d^{2}s_{2}/dz^{2}=\partial V/\partial s_{2}$, is
automatically satisfied when $s_{2}=0$. In this case, the domain wall solution
can be obtained by solving the equation of motion for $s_{1}$.

However, in the more general case without the alignment limit, it is not
necessary that only one Higgs boson gets VEV everywhere in the domain wall
solution for 2HDM. In that case, the gradient of the potential along the
$S_{1}$ direction ($s_{2}=0$) is

$\displaystyle\left.\frac{\partial V}{\partial s_{1}}\right|_{s_{2}=0}$
$\displaystyle=$
$\displaystyle\frac{m_{h}^{2}+m_{H}^{2}+(m_{H}^{2}-m_{h}^{2})\cos
2(\alpha-\beta)}{2v_{\rm SM}^{2}}s_{1}(s_{1}^{2}-v_{\rm SM}^{2})\,,$ (32)
$\displaystyle\left.\frac{\partial V}{\partial s_{2}}\right|_{s_{2}=0}$
$\displaystyle=$ $\displaystyle\frac{(m_{H}^{2}-m_{h}^{2})\sin
2(\alpha-\beta)}{2v_{\rm SM}^{2}}s_{1}(s_{1}^{2}-v_{\rm SM}^{2})\,.$ (33)

The equation of motion for $s_{2}$ is no longer trivial when $s_{2}=0$. In
order to get the correct domain wall solution, the path has to depart from the
$S_{1}$ direction. Although only $S_{1}$ gets a VEV at the minima, the value
of $S_{2}$ is non-zero inside the domain wall.

##  Appendix C Signal-to-noise ratio

We calculate the signal-to-noise ratio (SNR) of the gravitational wave signal
from the hybrid topological defect at various interferometers as follows [82].

SNR $\displaystyle=\sqrt{\tau\int_{f_{\rm min}}^{f_{\rm
max}}df\left[\frac{\Omega_{\rm GW}h^{2}(f)}{\Omega_{\rm
Sens}h^{2}}\right]^{2}},$ (34)

where $\Omega_{\rm Sens}h^{2}$ denotes the sensitivity of the interferometers,
$\Omega_{\rm GW}h^{2}$ is the GW signal, and $\tau$ is the duration of running
in years. We take $\tau=4$ for all interferometers. In this definition,
$\text{SNR}>10$ is taken to be the threshold for detection in a given
interferometer [82]. The resulting SNR is shown in Table 3.

Interferometer |  ${v_{\rm{M}}}=10^{12}$ GeV |  $v_{\rm{M}}=10^{13}$ GeV |  $v_{\rm{M}}=10^{14}$ GeV |  $v_{\rm{M}}=10^{15}$ GeV
---|---|---|---|---
$\mu$ARES | 0.0 | 0.0 | $6.9\times 10^{3}$ | $1.15\times 10^{6}$
LISA | 0.0 | 0.01 | 343.51 | $10^{4}$
BBO | 3.46 | $3.05\times 10^{5}$ | $1.28\times 10^{7}$ | $1.36\times 10^{8}$
DECIGO | 2.43 | $7.13\times 10^{4}$ | $3.51\times 10^{6}$ | $3.85\times 10^{7}$
AION | 0.0 | 0.97 | 57.51 | 630.99
AEDGE | 0.0 | 15.19 | $2.08\times 10^{3}$ | $2.61\times 10^{4}$
ET | 7.21 | 379.58 | $4.33\times 10^{3}$ | $3.78\times 10^{4}$
CE | 12.79 | 654.21 | $7.31\times 10^{3}$ | $6.38\times 10^{4}$
A+ | 0.01 | 0.24 | 2.58 | 22.55
HLV | 0.0 | 0.07 | 0.71 | 6.25
Table 3: Signal-to-noise ratio for $v_{\rm SM}=246$ GeV.

## References

  * [1] ATLAS collaboration, _Observation of a new particle in the search for the Standard Model Higgs boson with the ATLAS detector at the LHC_ , _Phys. Lett. B_ 716 (2012) 1 [1207.7214].
  * [2] CMS collaboration, _Observation of a New Boson at a Mass of 125 GeV with the CMS Experiment at the LHC_ , _Phys. Lett. B_ 716 (2012) 30 [1207.7235].
  * [3] CMS collaboration, _Search for a standard model-like Higgs boson in the mass range between 70 and 110 GeV in the diphoton final state in proton-proton collisions at $\sqrt{s}=$ 8 and 13 TeV_, _Phys. Lett. B_ 793 (2019) 320 [1811.08459].
  * [4] ILC collaboration, H. Baer et al., eds., _The International Linear Collider Technical Design Report - Volume 2: Physics_ , 1306.6352.
  * [5] F. An et al., _Precision Higgs physics at the CEPC_ , _Chin. Phys. C_ 43 (2019) 043002 [1810.09037].
  * [6] M. Cepeda et al., _Report from Working Group 2: Higgs Physics at the HL-LHC and HE-LHC_ , _CERN Yellow Rep. Monogr._ 7 (2019) 221 [1902.00134].
  * [7] T.D. Lee, _A Theory of Spontaneous T Violation_ , _Phys. Rev. D_ 8 (1973) 1226.
  * [8] S.L. Glashow and S. Weinberg, _Natural Conservation Laws for Neutral Currents_ , _Phys. Rev. D_ 15 (1977) 1958.
  * [9] E.A. Paschos, _Diagonal Neutral Currents_ , _Phys. Rev. D_ 15 (1977) 1966.
  * [10] G.C. Branco, P.M. Ferreira, L. Lavoura, M.N. Rebelo, M. Sher and J.P. Silva, _Theory and phenomenology of two-Higgs-doublet models_ , _Phys. Rept._ 516 (2012) 1 [1106.0034].
  * [11] I.P. Ivanov, _Building and testing models with extended Higgs sectors_ , _Prog. Part. Nucl. Phys._ 95 (2017) 160 [1702.03776].
  * [12] P. Ko, Y. Omura and C. Yu, _A Resolution of the Flavor Problem of Two Higgs Doublet Models with an Extra $U(1)_{H}$ Symmetry for Higgs Flavor_, _Phys. Lett. B_ 717 (2012) 202 [1204.4588].
  * [13] P. Ko, Y. Omura and C. Yu, _Higgs phenomenology in Type-I 2HDM with $U(1)_{H}$ Higgs gauge symmetry_, _JHEP_ 01 (2014) 016 [1309.7156].
  * [14] P. Ko, Y. Omura and C. Yu, _Dark matter and dark force in the type-I inert 2HDM with local U(1) H gauge symmetry_, _JHEP_ 11 (2014) 054 [1405.2138].
  * [15] M.D. Campos, D. Cogollo, M. Lindner, T. Melo, F.S. Queiroz and W. Rodejohann, _Neutrino Masses and Absence of Flavor Changing Interactions in the 2HDM from Gauge Principles_ , _JHEP_ 08 (2017) 092 [1705.05388].
  * [16] P. Minkowski, _$\mu\to e\gamma$ at a Rate of One Out of $10^{9}$ Muon Decays?_, _Phys. Lett. B_ 67 (1977) 421.
  * [17] T. Yanagida, _Horizontal gauge symmetry and masses of neutrinos_ , _$\mathrm{in}$  _Proceedings of the Workshop on Unified Theory and Baryon Number of the Universe_ , $\mathrm{Tsukuba,\ Japan}$_ (1979) .
  * [18] M. Gell-Mann, P. Ramond and R. Slansky, __in Sanibel talk_ , _retroprinted as_ _[9809459]_, _and in_ _Supergravity_ , _North-Holland, Amsterdam (1979), PRINT-80-0576, retroprinted as[1306.4669]__, .
  * [19] S. Glashow, _The Future of Elementary Particle Physics_ , _NATO Sci. Ser. B_ 61 (1980) 687.
  * [20] R.N. Mohapatra and G. Senjanovic, _Neutrino Mass and Spontaneous Parity Nonconservation_ , _Phys. Rev. Lett._ 44 (1980) 912.
  * [21] A. Arhrib, R. Benbrik and S. Moretti, _Bosonic Decays of Charged Higgs Bosons in a 2HDM Type-I_ , _Eur. Phys. J. C_ 77 (2017) 621 [1607.02402].
  * [22] A. Arhrib, R. Benbrik, S. Moretti, A. Rouchad, Q.-S. Yan and X. Zhang, _Multi-photon production in the Type-I 2HDM_ , _JHEP_ 07 (2018) 007 [1712.05332].
  * [23] R. Enberg, W. Klemm, S. Moretti and S. Munir, _Signatures of the Type-I 2HDM at the LHC_ , _PoS_ CORFU2018 (2018) 013 [1812.08623].
  * [24] Y. Wang, A. Arhrib, R. Benbrik, M. Krab, B. Manaut, S. Moretti et al., _Analysis of W ± \+ 4$\gamma$ in the 2HDM Type-I at the LHC_, _JHEP_ 12 (2021) 021 [2107.01451].
  * [25] Y. Wang, A. Arhrib, R. Benbrik, M. Krab, B. Manaut, S. Moretti et al., _Searching for $H^{\pm}\to W^{\pm}+4\gamma$ signals in the 2HDM Type-I at the LHC_, 11, 2021 [2111.12286].
  * [26] S. Moretti, S. Semlali and C.H. Shepherd-Themistocleous, _A novel experimental search channel for very light higgs bosons in the 2HDM type I_ , _Eur. Phys. J. C_ 84 (2024) 245 [2207.03007].
  * [27] S. Moretti, S. Semlali and C.H. Shepherd-Themistocleous, _Hunting light Higgses at the LHC in the context of the 2HDM Type-I_ , _PoS_ ICHEP2022 (2022) 529 [2211.11388].
  * [28] T. Mondal, S. Moretti, S. Munir and P. Sanyal, _Electroweak Multi-Higgs Production: A Smoking Gun for the Type-I Two-Higgs-Doublet Model_ , _Phys. Rev. Lett._ 131 (2023) 231801 [2304.07719].
  * [29] A. Arhrib, S. Moretti, S. Semlali, C.H. Shepherd-Themistocleous, Y. Wang and Q.S. Yan, _Searching for H→hh→bb¯ $\tau$$\tau$ in the 2HDM type-I at the LHC_, _Phys. Rev. D_ 109 (2024) 055020 [2310.02736].
  * [30] W. Liu, L. Wang and Y. Zhang, _Direct Production of Light Scalar in the Type-I Two-Higgs-Doublet Model at the Lifetime Frontier of LHC_ , 2403.16623.
  * [31] N. Chen, T. Han, S. Li, S. Su, W. Su and Y. Wu, _Type-I 2HDM under the Higgs and Electroweak Precision Measurements_ , _JHEP_ 08 (2020) 131 [1912.01431].
  * [32] T. Han, S. Li, S. Su, W. Su and Y. Wu, _Heavy Higgs bosons in 2HDM at a muon collider_ , _Phys. Rev. D_ 104 (2021) 055029 [2102.08386].
  * [33] T. Han, S. Li, S. Su, W. Su and Y. Wu, _Comparative Studies of 2HDMs under the Higgs Boson Precision Measurements_ , _JHEP_ 01 (2021) 045 [2008.05492].
  * [34] L. Bian, N. Chen, W. Su, Y. Wu and Y. Zhang, _Future prospects of mass-degenerate Higgs bosons in the CP -conserving two-Higgs-doublet model_ , _Phys. Rev. D_ 97 (2018) 115007 [1712.01299].
  * [35] N. Chen, T. Li, Z. Teng and Y. Wu, _Collapsing domain walls in the two-Higgs-doublet model and deep insights from the EDM_ , _JHEP_ 10 (2020) 081 [2006.06913].
  * [36] L. Delle Rose, S. Khalil and S. Moretti, _Explanation of the 17 MeV Atomki anomaly in a U(1)’ -extended two Higgs doublet model_ , _Phys. Rev. D_ 96 (2017) 115024 [1704.03436].
  * [37] ATLAS collaboration, _A detailed map of Higgs boson interactions by the ATLAS experiment ten years after the discovery_ , _Nature_ 607 (2022) 52 [2207.00092].
  * [38] CMS collaboration, _A portrait of the Higgs boson by the CMS experiment ten years after the discovery._ , _Nature_ 607 (2022) 60 [2207.00043].
  * [39] J.F. Gunion and H.E. Haber, _The CP conserving two Higgs doublet model: The Approach to the decoupling limit_ , _Phys. Rev. D_ 67 (2003) 075019 [hep-ph/0207010].
  * [40] H.E. Haber, _A natural mechanism for a SM-like Higgs boson in the 2HDM without decoupling_ , _PoS_ DISCRETE2020-2021 (2022) 010 [2205.07578].
  * [41] K. Saikawa, _A review of gravitational waves from cosmic domain walls_ , _Universe_ 3 (2017) 40 [1703.02576].
  * [42] C.L. Wainwright, _CosmoTransitions: Computing Cosmological Phase Transition Temperatures and Bubble Profiles with Multiple Fields_ , _Comput. Phys. Commun._ 183 (2012) 2006 [1109.4189].
  * [43] E.W. Kolb and M.S. Turner, eds., _THE EARLY UNIVERSE. REPRINTS_ (1988).
  * [44] R. Apreda, M. Maggiore, A. Nicolis and A. Riotto, _Gravitational waves from electroweak phase transitions_ , _Nucl. Phys. B_ 631 (2002) 342 [gr-qc/0107033].
  * [45] T.W.B. Kibble, _Topology of Cosmic Domains and Strings_ , _J. Phys. A_ 9 (1976) 1387.
  * [46] G.S.F. Guedes, P.P. Avelino and L. Sousa, _Signature of inflation in the stochastic gravitational wave background generated by cosmic string networks_ , _Phys. Rev. D_ 98 (2018) 123505 [1809.10802].
  * [47] Y. Cui, M. Lewicki and D.E. Morrissey, _Gravitational Wave Bursts as Harbingers of Cosmic Strings Diluted by Inflation_ , _Phys. Rev. Lett._ 125 (2020) 211302 [1912.08832].
  * [48] F. Ferrer, A. Ghoshal and M. Lewicki, _Imprints of a supercooled phase transition in the gravitational wave spectrum from a cosmic string network_ , _JHEP_ 09 (2023) 036 [2304.02636].
  * [49] T.W.B. Kibble, G. Lazarides and Q. Shafi, _Strings in SO(10)_ , _Phys. Lett. B_ 113 (1982) 237.
  * [50] T.W.B. Kibble, G. Lazarides and Q. Shafi, _Walls Bounded by Strings_ , _Phys. Rev. D_ 26 (1982) 435.
  * [51] L. Leblond, B. Shlaer and X. Siemens, _Gravitational Waves from Broken Cosmic Strings: The Bursts and the Beads_ , _Phys. Rev. D_ 79 (2009) 123519 [0903.4686].
  * [52] A.E. Everett and A. Vilenkin, _Left-right Symmetric Theories and Vacuum Domain Walls and Strings_ , _Nucl. Phys. B_ 207 (1982) 43.
  * [53] D.I. Dunsky, A. Ghoshal, H. Murayama, Y. Sakakihara and G. White, _GUTs, hybrid topological defects, and gravitational waves_ , _Phys. Rev. D_ 106 (2022) 075030 [2111.08750].
  * [54] X. Martin and A. Vilenkin, _Gravitational wave background from hybrid topological defects_ , _Phys. Rev. Lett._ 77 (1996) 2879 [astro-ph/9606022].
  * [55] J.J. Blanco-Pillado and K.D. Olum, _Stochastic gravitational wave background from smoothed cosmic string loops_ , _Phys. Rev. D_ 96 (2017) 104046 [1709.02693].
  * [56] J.J. Blanco-Pillado, K.D. Olum and B. Shlaer, _The number of cosmic string loops_ , _Phys. Rev. D_ 89 (2014) 023512 [1309.6637].
  * [57] P. Auclair et al., _Probing the gravitational wave background from cosmic strings with LISA_ , _JCAP_ 04 (2020) 034 [1909.00819].
  * [58] T. Vachaspati and A. Vilenkin, _Gravitational Radiation from Cosmic Strings_ , _Phys. Rev. D_ 31 (1985) 3052.
  * [59] G. Janssen et al., _Gravitational wave astronomy with the SKA_ , _PoS_ AASKA14 (2015) 037 [1501.00127].
  * [60] A. Sesana et al., _Unveiling the gravitational universe at $\mu$-Hz frequencies_, _Exper. Astron._ 51 (2021) 1333 [1908.11391].
  * [61] LISA collaboration, _Laser Interferometer Space Antenna_ , 1702.00786.
  * [62] H. Kudoh, A. Taruya, T. Hiramatsu and Y. Himemoto, _Detecting a gravitational-wave background with next-generation space interferometers_ , _Phys. Rev. D_ 73 (2006) 064006 [gr-qc/0511145].
  * [63] S. Kawamura et al., _Current status of space gravitational wave antenna DECIGO and B-DECIGO_ , _PTEP_ 2021 (2021) 05A105 [2006.13545].
  * [64] G.M. Harry, P. Fritschel, D.A. Shaddock, W. Folkner and E.S. Phinney, _Laser interferometry for the big bang observer_ , _Class. Quant. Grav._ 23 (2006) 4887.
  * [65] AEDGE collaboration, _AEDGE: Atomic Experiment for Dark Matter and Gravity Exploration in Space_ , _EPJ Quant. Technol._ 7 (2020) 6 [1908.00802].
  * [66] L. Badurina et al., _AION: An Atom Interferometer Observatory and Network_ , _JCAP_ 05 (2020) 011 [1911.11755].
  * [67] LIGO Scientific collaboration, _Exploring the Sensitivity of Next Generation Gravitational Wave Detectors_ , _Class. Quant. Grav._ 34 (2017) 044001 [1607.08697].
  * [68] S. Hild, S. Chelkowski and A. Freise, _Pushing towards the ET sensitivity using ’conventional’ technology_ , 0810.0604.
  * [69] LIGO Scientific, VIRGO, KAGRA collaboration, _Search for gravitational-wave transients associated with magnetar bursts in Advanced LIGO and Advanced Virgo data from the third observing run_ , 2210.10931.
  * [70] KAGRA, Virgo, LIGO Scientific collaboration, _Upper limits on the isotropic gravitational-wave background from Advanced LIGO and Advanced Virgo’s third observing run_ , _Phys. Rev. D_ 104 (2021) 022004 [2101.12130].
  * [71] Y. Jiang and Q.-G. Huang, _Upper limits on the polarized isotropic stochastic gravitational-wave background from advanced LIGO-Virgo’s first three observing runs_ , _JCAP_ 02 (2023) 026 [2210.09952].
  * [72] NANOGrav collaboration, _The NANOGrav 15 yr Data Set: Evidence for a Gravitational-wave Background_ , _Astrophys. J. Lett._ 951 (2023) L8 [2306.16213].
  * [73] NANOGrav collaboration, _The NANOGrav 15 yr Data Set: Search for Signals from New Physics_ , _Astrophys. J. Lett._ 951 (2023) L11 [2306.16219].
  * [74] EPTA collaboration, _The second data release from the European Pulsar Timing Array - I. The dataset and timing analysis_ , _Astron. Astrophys._ 678 (2023) A48 [2306.16224].
  * [75] EPTA, InPTA: collaboration, _The second data release from the European Pulsar Timing Array - III. Search for gravitational wave signals_ , _Astron. Astrophys._ 678 (2023) A50 [2306.16214].
  * [76] I. Baldes and M.O. Olea-Romacho, _Primordial black holes as dark matter: interferometric tests of phase transition origin_ , _JHEP_ 01 (2024) 133 [2307.11639].
  * [77] J.A. Dror, T. Hiramatsu, K. Kohri, H. Murayama and G. White, _Testing the Seesaw Mechanism and Leptogenesis with Gravitational Waves_ , _Phys. Rev. Lett._ 124 (2020) 041804 [1908.03227].
  * [78] A. Dasgupta, P.S.B. Dev, A. Ghoshal and A. Mazumdar, _Gravitational wave pathway to testable leptogenesis_ , _Phys. Rev. D_ 106 (2022) 075027 [2206.07032].
  * [79] B. Fu, A. Ghoshal and S.F. King, _Cosmic string gravitational waves from global $U(1)_{B-L}$ symmetry breaking as a probe of the type I seesaw scale_, _JHEP_ 11 (2023) 071 [2306.07334].
  * [80] P. Di Bari, S.F. King and M.H. Rahat, _Gravitational waves from phase transitions and cosmic strings in neutrino mass models with multiple Majorons_ , 2306.04680.
  * [81] S.F. King, D. Marfatia and M.H. Rahat, _Toward distinguishing Dirac from Majorana neutrino mass with gravitational waves_ , _Phys. Rev. D_ 109 (2024) 035014 [2306.05389].
  * [82] C. Caprini et al., _Science with the space-based interferometer eLISA. II: Gravitational waves from cosmological phase transitions_ , _JCAP_ 04 (2016) 001 [1512.06239].


# A matrix for counting paths in acyclic colored digraphs

Sudip Bera  Department of Mathematics, Indian Institute of Science, Bangalore
560 012<EMAIL_ADDRESS>
###### Abstract.

In this paper, we introduce a matrix $A_{\Gamma_{k,R}}$ associate with a
$k$-colored acyclic digraph $\Gamma_{k,R}$ such that
$\text{det}(A_{\Gamma_{k,R}})$ enumerates the paths in the digraph
$\Gamma_{k,R}.$

###### Key words and phrases:

Determinant; Colored digraph; Enumeration of paths

###### 2010 Mathematics Subject Classification:

05A19; 05A05; 05C30; 05C38

##  1\. Introduction

The evaluation of determinants is a nice topic, and fascinating for many
people [4, 7, 8, 9, 12]. Moreover, the evaluation of determinants in a
combinatorial way gives more insight into ‘why’ rather than ‘how’. This is
especially true when the entries of a matrix have natural graph theoretic
descriptions. Recently, in [1] the author shows a combinatorial interpretation
of the determinant of a matrix as a generating function over Brauer diagrams
in two different but related ways. This interpretation naturally explains why
the determinant of an even anti-symmetric matrix is the square of a Pfaffian.
In [2], the author gives a combinatorial interpretation of product formulas of
various types of hook determinants by counting weighted paths in a digraph.
Moreover in [11], Doron Zeilberger gives combinatorial proofs of some
determinantal identities. In [10], R. P. Stanley defines a matrix $A$
associated with an acyclic digraph $\Gamma$ (without multiple edges) such that
the determinant of $A$ enumerates the paths in $\Gamma.$ In fact, the theorem
of R. P. Stanley is as follows:

Let $[n]=\\{1,2,\cdots,n\\}$ be the vertex set of $\Gamma.$ Since $\Gamma$ is
acyclic, we may assume without loss of generality that $i<j$ whenever $(i,j)$
is an edge. Let $\\{x_{1},x_{2},\cdots,x_{n}\\}$ be a set of indeterminates.
Construct an $n\times n$ matrix $A=(a_{ij})$ as follows:

(1) $a_{ij}=\begin{cases}1+x_{i},&\text{ if }i=j\\\ x_{i},&\text{ if
}i<j\text{ and }(i,j)\text{ is not an edge }\\\ 0,&\text{ if }i<j\text{ and
}(i,j)\text{ is an edge }\\\ x_{i},&\text{ if }i>j.\end{cases}$

Then the theorem of Stanley says the following:

######  Theorem 1 ([10]).

${\rm det(A)}=1+\sum\limits_{P}x_{k_{1}}x_{k_{2}}\cdots x_{k_{j}}$, where $P$
ranges over all paths $k_{1}k_{2}\cdots k_{j}.$

The above theorem is an easy consequence of a more general theorem of Goulden
and Jackson [6]. Motivated by Theorem 1, in this article we define a matrix
$A_{\Gamma_{k,R}}$ associated with an acyclic digraph $\Gamma_{k,R}$ (having
multiple edges) such that $\text{det}(A_{\Gamma_{k,R}})$ enumerates the paths
in the digraph $\Gamma_{k,R}.$

##  2\. Basic definitions and main theorem

In this section, we state our main theorem. Before that, let us briefly
describe a well known graphical interpretation of determinant. See [5] for
details. To an $n\times n$ matrix $A,$ we can associate a weighted digraph
$D(A),$ with vertex set $[n]$ and for each ordered pair $(i,j),$ there is an
edge directed from $i$ to $j$ with weight $a_{ij}$. A _linear subdigraph_
$\gamma$ of $D(A)$ is a spanning collection of pairwise vertex-disjoint
cycles. A loop around a single vertex is also considered to be a cycle of
length $1.$ The weight of a linear subdigraph $\gamma,$ written as $w(\gamma)$
is the product of the weights of all its edges. The number of cycles contained
in $\gamma$ is denoted by $c(\gamma).$ Now the cycle-decomposition of
permutations yields the following description of ${\rm det}(A),$ namely

${\rm det}(A)=\sum\limits_{\gamma}(-1)^{n+c(\gamma)}w(\gamma),$

where the summation runs over all linear subdigraphs $\gamma$ of $D(A).$

The _standard representation_ of a directed cycle is defined to be a sequence
of vertices $v_{0}\rightarrow v_{1}\rightarrow\cdots\rightarrow
v_{k}\rightarrow v_{0}$ such that $v_{0}$ is the smallest among the elements
$\\{v_{0},v_{1},\cdots,v_{k}\\}.$ Although any vertex on a directed cycle can
be treated as the initial and terminal point of the cycle, in the standard
representation only the smallest vertex is treated to be the initial and the
terminal point of the cycle. A path $v_{i}\rightarrow v_{k_{1}}\rightarrow
v_{k_{2}}\rightarrow\cdots\rightarrow v_{k_{l}}\rightarrow v_{j}$ in a graph
$\Gamma$ is said to be _vertex-increasing_ if
$v_{i}<v_{k_{1}}<\cdots<v_{k_{\ell}}<v_{j}.$ Vertex-increasing path between
the vertices $v_{i},v_{j}$ is denoted by $v_{i}\text{I}v_{j}.$ Similarly we
may define a _vertex-decreasing_ path and denote the vertex-decreasing path
between the vertices $v_{i},v_{j}$ by $v_{i}\text{D}v_{j}.$ A path is said to
be _IDI_ path if the standard representation of this path contains a subpath
of the form $v_{i}\text{I}v_{j}\text{D}v_{k}\text{I}v_{\ell},$ for some
$v_{i},v_{j},v_{k},v_{\ell}.$ A vertex $v$ of a cycle $C_{i}$ in a linear
subdigraph is said to be _enclosed_ if it is surrounded by another cycle
$C_{j}(\neq C_{i})$ of that linear subdigraph. For example, in the first
diagram of Figure 2, $1\text{I}7\text{D}6\text{I}9$ is a IDI path and the
vertices $4,5$ are surrounded by the cycle $1\rightarrow 7\rightarrow
6\rightarrow 9\rightarrow 10\rightarrow 3\rightarrow 2\rightarrow 1.$ Hence
$4,5$ are enclosed vertices.

######  Definition 2.

A linear subdigraph is said to be _complex_ if and only if

  1. (1)

either at least one of its cycle encloses a vertex of another nontrivial
(i.e., not a loop) cycle

  2. (2)

or it has at least one cycle $C,$ whose standard representation contains a
subpath of the form $v_{i}\text{I}v_{j}\text{D}v_{k}\text{I}v_{\ell},$ (i.e.,
$C$ has a IDI path) for some $v_{i},v_{j},v_{k},v_{\ell}.$

A cycle $C$ of a complex linear subdigraph is said to be _singular_ if and
only either $C$ encloses a vertex of another nontrivial cycle of the same
linear subdigraph, or it has a IDI path (of the form
$v_{i}\text{I}v_{j}\text{D}v_{k}\text{I}v_{\ell},$ for some
$v_{i},v_{j},v_{k},v_{\ell}).$

######  Definition 3.

A vertex $v_{k}$ is said to be a _point of singularity_ of a singular cycle
$C$ if $v_{k}$ is a vertex of another nontrivial cycle of the same linear
subdigraph enclosed by $C$ or there exists a IDI path
$v_{i}\text{I}v_{j}\text{D}v_{k}\text{I}v_{\ell}$ along $C$ (for the latter
case, sometime we call $v_{k}$ a corner point of the IDI path).

For example, in Figure 2, the first diagram has three cycles namely,
$C_{1}=1\rightarrow 7\rightarrow 6\rightarrow 9\rightarrow 10\rightarrow
3\rightarrow 2\rightarrow 1,C_{2}=4\rightarrow 5\rightarrow 4\text{ and
}C_{3}=8\rightarrow 11\rightarrow 8.$ The point of singularity of $C_{1}$ are
$6$ (as 1I7D6I9 is a IDI path); $4,5$ (as $4,5$ are vertices of the cycle
$4\rightarrow 5\rightarrow 4\text{ and }4,5,\text{ are enclosed by }C_{1});$
and $8$ (as $8$ is a vertex of the cycle $8\rightarrow 11\rightarrow 8$ and
$8$ is enclosed by the cycle $C_{1}).$ Also, the point of singularity of
$C_{3}$ are $9,10$ (as $9,10$ are enclosed by the cycle $8\rightarrow
11\rightarrow 8).$ But $C_{2}$ has no point of singularity. So, here
$C_{1},C_{3}$ are singular cycles but $C_{2}$ is not a singular cycle.

######  Definition 4.

A _$k$ -colored digraph_, denoted by $\Gamma_{k,R}$ is a digraph equipped with
a finite set $R=\\{c_{1},c_{2},\cdots,c_{k}\\}$ called the set of _colors_
such that for any ordered pair of vertices $(i,j)$ in $\Gamma_{k,R},$ there
are $\ell(0\leq\ell\leq k)$ directed edges from $i$ to $j$ each receiving
distinct colors from $R$ (i.e., no two of the directed edges from $i$ to $j$
receive the same color).

A $k$-colored digraph $\Gamma_{k,R}$ is said to be _acyclic_ if it has no
cycle. In this case, without loss of generality we assume that $i<j$ whenever
there is an edge from $i$ to $j.$ i.e., in $\Gamma_{k,R}$ there is no edge
from the vertex $j$ to the vertex $i,$ for any $i<j.$ By
$v_{k_{1}}\underbrace{\rightarrow}_{c_{r_{1}}}v_{k_{2}}\underbrace{\rightarrow}_{c_{r_{2}}}v_{k_{3}}\underbrace{\rightarrow}_{c_{r_{3}}}v_{k_{4}}\underbrace{\rightarrow}_{c_{r_{4}}}\cdots\underbrace{\rightarrow}_{c_{r_{t-2}}}v_{k_{t-1}}\underbrace{\rightarrow}_{c_{r_{t-1}}}v_{k_{t}},$
we mean a color (not necessarily distinct) path from the vertex $v_{k_{1}}$ to
the vertex $v_{k_{t}}$ such that the edge from $v_{k_{i}}$ to $v_{k_{i+1}}$ is
colored by the color $c_{r_{i}},$ for $i=1,\cdots,t-1.$ Let the vertex set of
the graph $\Gamma_{k,R}$ be $[n].$ For each $i\in[n],$ we assign $k$ variables
$x_{i}^{(1)},x_{i}^{(2)},\cdots,x_{i}^{(k)},$ where for each
$j\in[k],x_{i}^{(j)}$ is the variable corresponding to the color $c_{j}.$ Let
$\\{x_{1}^{(1)},x_{1}^{(2)},\cdots,x_{1}^{(k)},x_{2}^{(1)},x_{2}^{(2)},\cdots,x_{2}^{(k)},\cdots,x_{n}^{(1)},x_{n}^{(2)},\cdots,x_{n}^{(k)}\\}$
be the set of variables. Let us consider an $n\times n$ matrix
$A_{\Gamma_{k,R}}=(a_{ij})$ as follows:

(2) $a_{ij}=\begin{cases}1+\sum\limits_{t\in[k]}x_{i}^{(t)},&\text{ if }i=j\\\
\sum\limits_{t\in[k]\setminus\\{i_{1},\cdots,i_{r}\\}}x_{i}^{(t)},&\text{ if
}i<j\text{ and }i\text{ to }j\text{ there are }r(0\leq r<k)\text{ different
}\\\ &\text{ colored edges using the colors }{c_{i_{1}}},\cdots,{c_{i_{r}}}\\\
0,&\text{ if }i<j\text{ and }i\text{ to }j\text{ there are }k\text{ different
colored edges }\\\ \sum\limits_{t\in[k]}x_{i}^{(t)},&\text{ if }i>j.\\\
\end{cases}$

Note that, if $i<j$ and $r=0,$ then the set
$\\{i_{1},\cdots,i_{r}\\}=\emptyset.$ So, in this case
$a_{ij}=\sum\limits_{t\in[k]}x_{i}^{(t)}.$ With the notation and terminology
introduced above we can state our main theorem.

######  Theorem 5.

Let $\Gamma_{k,R}$ be a $k$-colored acyclic digraph such that the vertex set
of $\Gamma_{k,R}$ is $[n].$ Let $A_{\Gamma_{k,R}}$ be an $n\times n$ matrix
defined as (2). Then $\text{det}(A_{\Gamma_{k,R}})=$

$\displaystyle
1+\sum\limits_{P}\left(x_{i_{1}}^{(r_{1})}x_{i_{2}}^{(r_{2})}\cdots
x_{i_{t-1}}^{(r_{t-1})}x_{i_{t}}^{(1)}+x_{i_{1}}^{(r_{1})}x_{i_{2}}^{(r_{2})}\cdots
x_{i_{t-1}}^{(r_{t-1})}x_{i_{t}}^{(2)}+\cdots+x_{i_{1}}^{(r_{1})}x_{i_{2}}^{(r_{2})}\cdots
x_{i_{t-1}}^{(r_{t-1})}x_{i_{t}}^{(k)}\right),$ $\displaystyle\text{ where
}P\text{ ranges over all paths
}i_{1}\underbrace{\rightarrow}_{c_{r_{1}}}i_{2}\underbrace{\rightarrow}_{c_{r_{2}}}i_{3}\underbrace{\rightarrow}_{c_{r_{3}}}\cdots\underbrace{\rightarrow}_{c_{r_{t-2}}}i_{t-1}\underbrace{\rightarrow}_{c_{r_{t-1}}}i_{t}\text{
in the graph }\Gamma_{k,R}.$

######  Remark 6.

Note that, if the color set $R$ contains exactly one color (that is if
$R=\\{c\\},$ where $c$ represents a color), then $\Gamma_{1,R}$ is acyclic
digraph without multiple edges. Therefore as a corollary of this theorem, we
get Theorem 1.

In this portion, we explain our main theorem with one example.

######  Example 7.

Here we first fix the color set $R=\\{c_{1},c_{2},c_{3}\\},$ where
$c_{1},c_{2},c_{3}$ represent the colors black, blue, red respectively. Also
we choose the variables
$x_{1}^{(1)},x_{1}^{(2)},x_{1}^{(3)},x_{2}^{(1)},x_{2}^{(2)},x_{2}^{(3)},x_{3}^{(1)},x_{3}^{(2)},x_{3}^{(3)},x_{4}^{(1)},x_{4}^{(2)},x_{4}^{(3)},$
where we use the variable $x_{i}^{(j)}$ for the edge coming out from the
vertex $i$ and colored by $c_{j}$-th color, $j\in[3].$

$\bf{v_{2}}$$\bf{v_{1}}$$\bf{v_{3}}$$\bf{v_{4}}$ Figure 1. A colored digraph
$\Gamma_{3,R}$

Now from Definition (2), the matrix corresponding to the graph $\Gamma$ in
Figure 1 is the following:

$A_{\Gamma_{3,R}}=\left(\begin{array}[]{cccc}1+\sum\limits_{i=1}^{3}x_{1}^{(i)}&x_{1}^{(3)}&\sum\limits_{i=1}^{3}x_{1}^{(i)}&\sum\limits_{i=1}^{3}x_{1}^{(i)}\\\
\sum\limits_{i=1}^{3}x_{2}^{(i)}&1+\sum\limits_{i=1}^{3}x_{2}^{(i)}&\sum\limits_{i=2}^{3}x_{2}^{(i)}&\sum\limits_{i=2}^{3}x_{2}^{(i)}\\\
\sum\limits_{i=1}^{3}x_{3}^{(i)}&\sum\limits_{i=1}^{3}x_{3}^{(i)}&1+\sum\limits_{i=1}^{3}x_{3}^{(i)}&\sum\limits_{i=1}^{3}x_{3}^{(i)}\\\
\sum\limits_{i=1}^{3}x_{4}^{(i)}&\sum\limits_{i=1}^{3}x_{4}^{(i)}&\sum\limits_{i=1}^{3}x_{4}^{(i)}&1+\sum\limits_{i=1}^{3}x_{4}^{(i)}\end{array}\right).$

$\text{ Then, }\text{det}(A_{\Gamma_{3,R}})$ is equal to

$\displaystyle 1+$ $\displaystyle
x_{1}^{(1)}+x_{1}^{(2)}+x_{1}^{(3)}+x_{2}^{(1)}+x_{2}^{(2)}+x_{2}^{(3)}+x_{3}^{(1)}+x_{3}^{(2)}+x_{3}^{(3)}+x_{4}^{(1)}+x_{4}^{(2)}+x_{4}^{(3)}$
$\displaystyle+$ $\displaystyle
x_{1}^{(1)}x_{2}^{(1)}+x_{1}^{(1)}x_{2}^{(2)}+x_{1}^{(1)}x_{2}^{(3)}+x_{1}^{(2)}x_{2}^{(1)}+x_{1}^{(2)}x_{2}^{(2)}+x_{1}^{(2)}x_{2}^{(3)}+x_{2}^{(1)}x_{3}^{(1)}+x_{2}^{(1)}x_{3}^{(2)}+x_{2}^{(1)}x_{3}^{(3)}$
$\displaystyle+$ $\displaystyle
x_{2}^{(1)}x_{4}^{(1)}+x_{2}^{(1)}x_{4}^{(2)}+x_{2}^{(1)}x_{4}^{(3)}+x_{1}^{(1)}x_{2}^{(1)}x_{3}^{(1)}+x_{1}^{(1)}x_{2}^{(1)}x_{3}^{(2)}+x_{1}^{(1)}x_{2}^{(1)}x_{3}^{(3)}+x_{1}^{(2)}x_{2}^{(1)}x_{3}^{(1)}$
$\displaystyle+$ $\displaystyle
x_{1}^{(2)}x_{2}^{(1)}x_{3}^{(2)}+x_{1}^{(2)}x_{2}^{(1)}x_{3}^{(3)}+x_{1}^{(1)}x_{2}^{(1)}x_{4}^{(1)}+x_{1}^{(1)}x_{2}^{(1)}x_{4}^{(2)}+x_{1}^{(1)}x_{2}^{(1)}x_{4}^{(3)}$
$\displaystyle+$ $\displaystyle
x_{1}^{(2)}x_{2}^{(1)}x_{4}^{(1)}+x_{1}^{(2)}x_{2}^{(1)}x_{4}^{(2)}+x_{1}^{(2)}x_{2}^{(1)}x_{4}^{(3)}.$

##  3\. Proof of main theorem

We now prove our main theorem in a combinatorial way. Here we use the idea of
the proof of Theorem 1 described in [3].

###### Proof of Theorem 5.

Take

$W=\\{x_{1}^{(1)},x_{1}^{(2)},\cdots,x_{1}^{(k)},x_{2}^{(1)},x_{2}^{(2)},\cdots,x_{2}^{(k)},\cdots,x_{n}^{(1)},x_{n}^{(2)},\cdots,x_{n}^{(k)}\\}$

as the set of _letters_. The _free monoid_ $W^{*}$ is the set of all finite
sequences (including the empty sequence, denoted by $1$) of elements of $W,$
usually called _words_ , with the operation of concatenation. Construct an
algebra from $W^{*}$ by taking formal sum of elements of $W$ with coefficient
in $\mathbb{Z},$ extending the multiplication by usual distributivity. For
example, in this algebra,

$\displaystyle(x_{i}^{(1)}+x_{j}^{(2)})(x_{i}^{(1)}+x_{j}^{(2)})=$
$\displaystyle(x_{i}^{(1)}+x_{j}^{(2)})^{2}=x_{i}^{(1)}x_{i}^{(1)}+x_{i}^{(1)}x_{j}^{(2)}+x_{j}^{(2)}x_{i}^{(1)}+x_{j}^{(2)}x_{j}^{(2)},$
$\displaystyle(1+x_{i}^{(1)})x_{j}^{(2)}=x_{j}^{(2)}+x_{i}^{(1)}x_{j}^{(2)},\text{
etc}.$

A nonempty word $x_{k_{1}}^{(i_{1})}x_{k_{2}}^{(i_{2})}\cdots
x_{k_{m}}^{(i_{m})}$ is said to be _nice_ if $k_{1}<k_{2}<\cdots<k_{m}.$ The
empty word $1$ is also considered to be nice. Now let $B$ be a set of ordered
pairs $(x_{k_{r}}^{(i_{r})},x_{k_{s}}^{(i_{s})})$ with $k_{r}<k_{s}.$ Each
pair in the set $B$ is called a _bad pair_. In terms of digraph,
$(x_{k_{r}}^{(i_{r})},x_{k_{s}}^{(i_{s})})$ is a bad pair if and only if there
is no $i_{r}$-color edge (an edge colored by $i_{r}$th-color) from the vertex
$k_{r}$ to the vertex $k_{s}.$ As a result,
$(x_{k_{r}}^{(i_{r})},x_{k_{s}}^{(i_{s})})$ is being a bad pair implies that

$(x_{k_{r}}^{(i_{r})},x_{k_{s}}^{(1)}),(x_{k_{r}}^{(i_{r})},x_{k_{s}}^{(2)}),\cdots,(x_{k_{r}}^{(i_{r})},x_{k_{s}}^{(k)})$

are also bad pairs. We call a word $w\in W^{*}$ _best_ if $w$ is nice and no
bad pair can occur in $w$ as a consecutive pair. Let $S$ be the sum of all
best words in $W^{*}.$ Note that $1$ is also a best word. Now we derive an
expression for $S$ by using the _Principle of Inclusion Exclusion_(in short
_PIE_). The sum of all possible nice words is

(3) $\displaystyle
W_{0}=\prod\limits_{g=1}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

The sum of all possible nice words, where there is an occurrence of at least
one bad pair $(x_{r}^{(i)},x_{s}^{(j)})$ (let $s=r+c,c>0$) as a consecutive
pair is

(4) $\displaystyle
w^{(i,j)}_{rs}=\prod\limits_{g=1}^{r-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)x_{r}^{(i)}x_{s}^{(j)}\prod\limits_{g=s+1}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

Let $W_{1}$ be the sum of all possible nice words, where there is an
occurrence of at least one bad pair as a consecutive pair. The sum of all
possible nice words, where there is an occurrence of at least two bad pairs
$(x_{p}^{(i)},x_{q}^{(j)})$ and $(x_{r}^{(\ell)},x_{s}^{(m)})$ (where
$(p<q,r<s)$ and $i,j,\ell,m\in[k]$) as two consecutive pairs is

(5) $\displaystyle
w^{(i,j,\ell,m)}_{pq,rs}=\prod\limits_{g=1}^{p-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)x_{p}^{(i)}x_{q}^{(j)}\prod\limits_{g=q+1}^{r-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)x_{r}^{(\ell)}x_{s}^{(m)}\prod\limits_{g=s+1}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

Let $W_{2}$ be the sum of all possible nice words, where there is an
occurrence of at least two bad pairs as consecutive pairs. Define
$W_{3},W_{4},\cdots$ etc., in a similar fashion. Then by the PIE,

$S=W_{0}-W_{1}+W_{2}-W_{3}+\cdots.$

Now, to any term of an arbitrary $W_{i}$ in the above expression, we associate
a term or terms of $\text{det}(A_{\Gamma_{k,R}})$ so that the weights and
signs match. Actually we want to show that the image of a word under the ring
homomorphism
$f:W^{*}\rightarrow\mathbb{Z}\left[x_{1}^{(1)},\cdots,x_{1}^{(k)},\cdots,x_{n}^{(1)},\cdots,x_{n}^{(k)}\right],f(x_{i}^{(j)})=x_{i}^{(j)}$
(for all $i,j)$ is the sum of the weights of all linear subdigraphs associated
to that word. This can be accomplished in the following way.

For the sake of simplicity, here we demonstrate, how to associate to each term
of $W_{0}$ and $W_{1}$ in the above expression, a term or terms of
$\text{det}(A_{\Gamma_{k,R}}).$ The case of other indices
$W_{2},W_{3},W_{4},\cdots$ are similar. First of all $W_{0}$ consists the
single term

(6)
$\displaystyle\prod\limits_{g=1}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

To this term, we associate the linear subdigraph consisting of $n$ disjoint
loops around each vertex. Note that

$a_{ii}=\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right),\text{ for each
}i\in[n].$

Consequently, the sign weight of this digraph is same as (6). Now let us
assume that in the digraph $\Gamma,$ there are exactly $m$-colored edges from
the vertex $r$ to the vertex $s$ receiving $i_{1},\cdots,i_{m}$ different
colors. Hence from the vertex $r$ to the vertex $s,$ there is no $i$-colored
edge, where $i\in[k]\setminus\\{i_{1}\cdots,i_{m}\\}.$ As a result, in the
matrix $A_{\Gamma_{k,C}}=(a_{i,j})_{n\times n},$
$a_{r,s}=\sum\limits_{\ell=m+1}^{k}x_{r}^{(i_{\ell})}.$ Consider the set

(7) $\displaystyle
W_{1}(r,s)=\\{\prod\limits_{i=1}^{r-1}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right)x_{r}^{(i_{\ell})}x_{s}^{(d)}\prod\limits_{i=s+1}^{n}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right):\ell\in[k]\setminus[m],d\in[k]\\}.$

Clearly, each member of the set $W_{1}(r,s)$ is a typical term in $W_{1}.$ The
sign weighted sum of all elements in the set ${W_{1}(r,s)}$ (defined as (7))is

(8)
$-\prod\limits_{i=1}^{r-1}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right)\sum\limits_{\ell=m+1}^{k}\sum\limits_{d=1}^{k}x_{r}^{(i_{\ell})}x_{s}^{(d)}\prod\limits_{i=s+1}^{n}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right)$

Now to (8), we associate $2^{(c-1)}$ linear subdigraphs such that the sum of
sign weights of $2^{(c-1)}$ linear subdigraphs is (8). First consider the set
$T=\\{r+1,r+2,\cdots,r+c-1\\}\subset[n].$ Let
$T^{\prime}=\\{t_{1},\cdots,t_{m}\\}$ be an arbitrary subset of $T$ such that
$t_{1}<t_{2}<\cdots<t_{m}.$ Here we associate a linear subdigraph
$L_{T^{\prime}}$ as follows;

$L_{T^{\prime}}$ consists of a cycle $r\rightarrow s\rightarrow
t_{m}\rightarrow t_{m-1}\rightarrow\cdots\rightarrow t_{2}\rightarrow
t_{1}\rightarrow r$ and $(n-m-2)$ disjoint loops around the remaining
$(n-m-2)$ vertices. So, for each non empty subset $T^{\prime},$ we get a
unique linear subdigraph $L_{T^{\prime}}$ and

$\displaystyle
w(L_{T^{\prime}})=(-1)^{m+1}\left(a_{rs}a_{st_{m}}a_{t_{m}t_{m-1}}\cdots
a_{t_{2}t_{1}}a_{t_{1}r}\right)\prod\limits_{i\in[n]\setminus T}a_{ii},$
$\displaystyle=(-1)^{m+1}$
$\displaystyle\left(\sum\limits_{\ell=m+1}^{k}x_{r}^{(i_{\ell})}\sum\limits_{a=1}^{k}x_{s}^{(a)}\sum\limits_{a=1}^{k}x_{t_{m}}^{(a)}\cdots\sum\limits_{a=1}^{k}x_{t_{2}}^{(a)}\sum\limits_{a=1}^{k}x_{t_{1}}^{(a)}\right)\prod\limits_{i\in[n]\setminus
T}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right).$

If $T^{\prime}=\varphi,$ then $L_{T^{\prime}}$ consists of a cycle
$r\rightarrow s\rightarrow r$ and $(n-2)$ disjoint loops around the remaining
$(n-2)$ vertices and

$\displaystyle
w(L_{\varphi})=-a_{rs}a_{sr}\prod\limits_{i\in[n]\setminus\\{r,s\\}}a_{ii}$
$\displaystyle=-$
$\displaystyle\sum\limits_{\ell=m+1}^{k}x_{r}^{(i_{\ell})}\sum\limits_{a=1}^{k}x_{s}^{(a)}\prod\limits_{i\in[n]\setminus\\{r,s\\}}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right)$

Now the sum of the sign weights of all these $2^{(c-1)}$ linear subdigraphs is

(9)
$-\prod\limits_{i=1}^{r-1}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right)\left(\sum\limits_{\ell=m+1}^{k}\sum\limits_{d=1}^{k}x_{r}^{(i_{\ell})}x_{s}^{(d)}\right)X\prod\limits_{i=s+1}^{n}\left(1+\sum\limits_{a=1}^{k}x_{i}^{(a)}\right),$

where

$\displaystyle
X=\prod\limits_{g=r+1}^{s-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)-\sum\limits_{i\in
T}\prod\limits_{g=r+1}^{i-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)\left(\sum\limits_{a=1}^{k}x_{i}^{(a)}\right)\prod\limits_{g=i+1}^{s-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)$
$\displaystyle+\sum\limits_{i<j\in
T}\prod\limits_{g=1}^{i-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)\left(\sum\limits_{a=1}^{k}x_{i}^{(a)}\right)\prod\limits_{g=i+1}^{j-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)\left(\sum\limits_{a=1}^{k}x_{j}^{(a)}\right)\prod\limits_{g=j+1}^{s-1}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)$
$\displaystyle-\cdots+(-1)^{s-r-1}\prod\limits_{g=r+1}^{s-1}\left(\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

Let
$\mathcal{W}=\\{x_{r+1}^{(1)},\cdots,x_{r+1}^{(k)},x_{r+2}^{(1)},\cdots,x_{r+2}^{(k)},\cdots,x_{s-1}^{(1)},\cdots,x_{s-1}^{(k)}\\}$
be the set of letters. Since the above expression is the PIE model for the sum
of all nice words in $\mathcal{W}^{*}$ consisting of no letters, $X$ is
precisely $1.$

For example, if $r=1$ and $s=2,$ then $W_{1}(1,2)$ contains the following
terms

$\displaystyle x_{1}^{(i_{m+1})}x_{2}^{(1)}$
$\displaystyle\prod\limits_{g=3}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right),\cdots,x_{1}^{(i_{m+1})}x_{2}^{(k)}\prod\limits_{g=3}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right),\cdots,$
$\displaystyle
x_{1}^{(i_{k})}x_{2}^{(1)}\prod\limits_{g=3}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right),\cdots,x_{1}^{(i_{k})}x_{2}^{(k)}\prod\limits_{g=3}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

It is easy to see that, the sign of each term of $W_{1}(1,2)$ in $S$ is $-.$
Moreover, all the members of the set $W_{1}(1,2)$ are some typical terms of
$W_{1}.$ Now the sign sum of all these terms (terms of $W_{1}(1,2)$) is

$-\sum\limits_{\ell=m+1}^{k}\sum\limits_{d=1}^{k}x_{1}^{(i_{\ell})}x_{2}^{(d)}\prod\limits_{g=3}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

In this case, we associate the linear subdigraph, consisting of the cycle
$1\rightarrow 2\rightarrow 1$ and $(n-2)$ disjoint loops around the remaining
$(n-2)$ vertices. So, the number of cycle in this linear subdigraph is $n-1.$
Therefore, the weight of this linear subdigraph is

$\sum\limits_{\ell=m+1}^{k}\sum\limits_{d=1}^{k}x_{1}^{(i_{\ell})}x_{2}^{(d)}\prod\limits_{g=3}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)$

(as in the matrix
$a_{1,2}=\sum\limits_{\ell=m+1}^{k}x_{1}^{(i_{\ell})},a_{2,1}=\sum\limits_{d=1}^{k}x_{2}^{(d)},\text{
and for }i\geq 3,a_{i,i}=1+\sum\limits_{d=1}^{k}x_{i}^{(d)})$ and its sign in
$\text{det}(A_{\Gamma_{k,R}})$ is $-.$

Now if $r=1$ and $s=4,$ then

$W_{1}(1,4)=\\{x_{1}^{(i_{\ell})}x_{4}^{(d)}\prod\limits_{g=5}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right):\ell\in\\{m+1,\cdots,k\\},d\in[k]\\}.$

Here each term of the set $W_{1}(1,4)$ is a typical term of $W_{1}.$ Also the
sign sum of all the terms in the set $W_{1}(1,4)$ is

(10)
$-\sum\limits_{\ell=m+1}^{k}\sum\limits_{d=1}^{k}x_{1}^{(i_{\ell})}x_{4}^{(d)}\prod\limits_{g=5}^{n}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right).$

Now to (10) we associate four linear subdigraphs such that the sum of the
signed weights of the four linear subdigraphs is same as (10).

The first one consists of a cycle $1\rightarrow 4\rightarrow 1$ and $(n-2)$
disjoint loops around the remaining $(n-2)$ vertices. The second one consists
of a cycle $1\rightarrow 4\rightarrow 2\rightarrow 1$ and $(n-3)$ disjoint
loops around the remaining $(n-3)$ vertices. The third one consists of a cycle
$1\rightarrow 4\rightarrow 3\rightarrow 1$ and $(n-3)$ disjoint loops around
the remaining $(n-3)$ vertices. The last one consists of a cycle $1\rightarrow
4\rightarrow 3\rightarrow 2\rightarrow 1$ and $(n-4)$ disjoint loops around
the remaining $(n-4)$ vertices. Note that the signed weights of these four
linear subdigraphs are

$\displaystyle-\sum\limits_{\ell=m+1}^{k}\sum\limits_{d^{\prime}=1}^{k}x_{1}^{(i_{\ell})}x_{4}^{(d^{\prime})}\prod\limits_{g\in[k]\setminus\\{1,4\\}}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right),$
$\displaystyle+\sum\limits_{\ell=m+1}^{k}\sum\limits_{d^{\prime}=1}^{k}\sum\limits_{c^{\prime}=1}^{k}x_{1}^{(i_{\ell})}x_{4}^{(d^{\prime})}x_{2}^{(c^{\prime})}\prod\limits_{g\in[k]\setminus\\{1,2,4\\}}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right),$
$\displaystyle+\sum\limits_{\ell=m+1}^{k}\sum\limits_{d^{\prime}=1}^{k}\sum\limits_{c=1}^{k}x_{1}^{(i_{\ell})}x_{4}^{(d^{\prime})}x_{3}^{(c^{\prime})}\prod\limits_{g\in[k]\setminus\\{1,3,4\\}}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right),$
$\displaystyle-\sum\limits_{\ell=m+1}^{k}\sum\limits_{d^{\prime}=1}^{k}\sum\limits_{c^{\prime}=1}^{k}\sum\limits_{b^{\prime}=1}^{k}x_{1}^{(i_{\ell})}x_{4}^{(d^{\prime})}x_{3}^{(c^{\prime})}x_{2}^{(b^{\prime})}\prod\limits_{g\in[k]\setminus\\{1,2,3,4\\}}\left(1+\sum\limits_{a=1}^{k}x_{g}^{(a)}\right)$

respectively. The contribution of these four linear subdigraphs to
$\text{det}(A_{\Gamma_{k,R}})$ is the same as (10), since

$\displaystyle\left(1+\sum\limits_{a=1}^{k}x_{2}^{(a)}\right)$
$\displaystyle\left(1+\sum\limits_{a=1}^{k}x_{3}^{(a)}\right)-\left(\sum\limits_{a=1}^{k}x_{2}^{(a)}\right)\left(1+\sum\limits_{a=1}^{k}x_{3}^{(a)}\right)-$
$\displaystyle\left(1+\sum\limits_{a=1}^{k}x_{2}^{(a)}\right)\left(\sum\limits_{a=1}^{k}x_{3}^{(a)}\right)+\left(\sum\limits_{a=1}^{k}x_{2}^{(a)}\right)\left(\sum\limits_{a=1}^{k}x_{3}^{(a)}\right)=1.$

So, we have associated some terms of $\text{det}(A_{\Gamma_{k,R}})$ to the
terms of $S,$ matching signs and weights. But in
$\text{det}(A_{\Gamma_{k,R}}),$ beside these linear subdigraphs, that have
been used so far, there are other linear subdigraphs, which are complex linear
subdigraphs according to Definition 2. Therefore, to prove the theorem, we
have to prove that the contribution of the complex linear subdigraphs to
$\text{det}(A_{\Gamma_{k,R}})$ is zero. We prove this by defining a weight
preserving and sign-reversing involution on $\mathcal{C},$ where $\mathcal{C}$
is the set of all complex linear subdigraphs of the digraph
$D(A_{\Gamma_{k,R}}).$ We now describe the involution.

First, consider a complex linear subdigraph. Choose the singular cycle in the
linear subdigraph with the smallest initial point (in the standard
representation). Find its smallest point of singularity. Either it is an
enclosed point or a corner point of a IDI path (by Definitions 2, 3). Now let
us define a map $f:\mathcal{C}\rightarrow\mathcal{C}$ such that $f$ is a sign
reversing and weight preserving involution. According to the definition of
complex linear subdigraph, we divide the definition of $f$ into two cases:

Case 1: Let $\gamma\in\mathcal{C}.$ Let the smallest point of singularity of
$\gamma$ be a enclosed point. Let the singular cycle with smallest initial
point be $C$ with the initial point $v_{1}.$ Suppose $v^{\prime}_{t}$ is the
smallest point of singularity of $C.$ i.e., $v^{\prime}_{t}$ is enclosed by
$C.$ Let $C^{\prime}=v^{\prime}_{t}\rightarrow\cdots\rightarrow
v^{\prime}_{r}\rightarrow v^{\prime}_{t}$ be the cycle with the initial point
$v^{\prime}_{t}.$ Since $v^{\prime}_{t}$ is enclosed by the cycle $C,$ then
there exists two consecutive vertices $v_{i},v_{j}$ in the cycle $C$ such that
$v_{i}<v^{\prime}_{t}<v_{j}$ and $v_{j}\rightarrow v_{i}$ (as $v_{t}^{\prime}$
is enclosed by the cycle $C).$ Now for this $\gamma,$ we define $f(\gamma)$ by
the following algorithm:

  1. (1)

First, we follow the path from $v_{1}$ to $v_{j}$ along the enclosing cycle
$C.$

  2. (2)

Then move from $v_{j}$ to $v^{\prime}_{t}$ by a directed edge
$v_{j}\rightarrow v^{\prime}_{t}$ (it is possible as $v^{\prime}_{t}<v_{j},$
so in $\Gamma_{k,R}$ there is no directed edge (as $\Gamma_{k,R}$ is acyclic)
from the vertex $v_{j}$ to the vertex $v^{\prime}_{t}.$ As a result, in
$A(\Gamma_{k,R}),a_{v_{j}v^{\prime}_{t}}\neq 0$ (by (2)) and hence in
$D(A(\Gamma_{k,R})),$ there is a an edge from $v_{j}$ to $v^{\prime}_{t}$ such
that weight of the edge is nonzero).

  3. (3)

Then follow the path from $v^{\prime}_{t}$ to $v^{\prime}_{r}$ along the cycle
$C^{\prime},$ where $v^{\prime}_{r}$ is a vertex in $C^{\prime}$ such that
$v^{\prime}_{t}<v^{\prime}_{r}$ and $v^{\prime}_{r}\rightarrow
v^{\prime}_{t}.$

  4. (4)

After that we move from the vertex $v^{\prime}_{r}$ to the vertex $v_{i}$ by a
directed edge $v^{\prime}_{r}\rightarrow v_{i}.$

  5. (5)

Finally, we shift from the vertex $v_{i}$ to the vertex $v_{1},$ along the
enclosing cycle $C.$

  6. (6)

All other cycles in $\gamma$ remain unchanged.

Therefore, we have a new cycle
$C^{\prime\prime}=v_{1}\rightarrow\cdots\rightarrow v_{j}\rightarrow
v^{\prime}_{t}\rightarrow\cdots\rightarrow v^{\prime}_{r}\rightarrow
v_{i}\rightarrow\cdots\rightarrow v_{1}$ by connecting the cycles $C$ and
$C^{\prime}.$ So, $f(\gamma)$ contains the cycle $C^{\prime\prime}$ and all
other cycles (in $\gamma$) which are neither $C$ nor $C^{\prime}.$ Clearly,
$f(\gamma)$ is a linear subdigraph. Also, it is possible to find a vertex
$v_{q}$ in $C$ such that there exist a longest vertex-decreasing path along
$C$ (possibly $v_{q}=v_{j})$ from the vertex $v_{q}$ to the vertex $v_{j}.$
Again, we can choose a vertex $v_{p}$ in $C$ such that $v_{p}$ to $v_{q}$ we
have a vertex-increasing path along $C$ (as $v_{q}$ and $v_{j}$ are connected
by a longest vertex-decreasing path and $v_{1}$ is the initial point of $C).$
Therefore, we get a vertex-increasing path from $v_{p}$ to $v_{q}$ and a
vertex-decreasing path from $v_{q}$ to $v^{\prime}_{t}$ (as
$v^{\prime}_{t}<v_{j}).$ Now $v^{\prime}_{t}$ is the smallest vertex among the
vertices of the cycle $C^{\prime}.$ So by the $3^{\text{rd}}$ part of the
definition of $f,$ we get vertex-increasing path from the vertex
$v^{\prime}_{t}$ to the vertex $v^{\prime}_{s}$ (say) along $C^{\prime}.$ As a
result, we get
$v_{p}\text{I}v_{q}\text{D}v^{\prime}_{t}\text{I}v^{\prime}_{s}$ in
$f(\gamma).$ So $v^{\prime}_{t}$ is a corner of a IDI path and
$v^{\prime}_{t}$ is the smallest singular point of the cycle
$C^{\prime\prime}.$ Hence $f(\gamma)\in\mathcal{C}.$ Moreover,
$w(\gamma)=w(f(\gamma))$ and in this case, the number of cycles in $f(\gamma)$
is decreased by $1.$

Case 2: Let the smallest point of singularity be a corner of a IDI path. Let
$\tilde{\gamma}\in\mathcal{C}$ such that $\tilde{C}$ is the singular cycle
with smallest initial point $u_{1}.$ Let $u^{\prime}_{t}$ be the singular
point such that $u^{\prime}_{t}$ be a corner of a IDI path
$\mathcal{P}=u_{1}\text{I}u_{j}\text{D}u^{\prime}_{t}\text{I}u_{r}.$ Now the
path $P$ from the vertex $u_{r}$ to $u_{1}$ has a maximal decreasing subpath
$P^{\prime}$ from a vertex $u_{s}$ (say) to $u_{1}.$ Clearly, on $P$ there are
two consecutive vertices $u_{p},u_{q}$ such that $u_{p}<u^{\prime}_{t}<u_{q}$
and $u_{q}\rightarrow u_{p}.$ Let $P^{\prime\prime}$ be a path starting from
$u_{j}$ such that $P^{\prime\prime}$ is same as the path
$u_{j}\text{D}u^{\prime}_{t}$ except the ending points and $P^{\prime\prime}$
does not contain the vertex $u^{\prime}_{t}.$ In fact, the ending points of
$u_{j}\text{D}u^{\prime}_{t},P^{\prime\prime}$ are $u^{\prime}_{t},u_{p}$
respectively. Now we follow the algorithm below for $f(\tilde{\gamma}).$

  1. (1)

First, follow the vertex-increasing path $u_{1}\text{I}u_{j}$ from $u_{1}$ to
$u_{j}$ along $\mathcal{P}.$

  2. (2)

Move from $u_{j}$ to $u_{p}$ along the path $P^{\prime\prime}$ and then shift
from $u_{p}$ to $u_{1}$ along $P^{\prime}$ to complete the cycle
$u_{1}\rightarrow\cdots\rightarrow u_{j}\rightarrow\cdots\rightarrow
u_{p}\rightarrow\cdots\rightarrow u_{1}.$

  3. (3)

After that from $u_{t}^{\prime}$ to $u_{r}$ move along the vertex-increasing
path $u^{\prime}_{t}\text{I}u_{r}$ on $\mathcal{P}.$

  4. (4)

Now follow the path $P$ from the vertex $u_{r}$ to the vertex $u_{q}.$

  5. (5)

Finally, move by a directed edge from $u_{q}$ to $u^{\prime}_{t}$ (it is
possible as $u^{\prime}_{t}<u_{q}$ and in $\Gamma$ there is no edge from
$u_{q}$ to $u^{\prime}_{t}$ and hence in
$A(\Gamma_{k,R}),a_{u_{q}u^{\prime}_{t}}\neq 0)$ to complete the cycle
$u^{\prime}_{t}\rightarrow\cdots\rightarrow u_{r}\rightarrow\cdots\rightarrow
u_{q}\rightarrow u^{\prime}_{t}.$

  6. (6)

All other cycles in $\tilde{\gamma}$ remain unchanged.

So, by our algorithm we split the cycle $\tilde{C}$ into two new vertex
disjoint cycles $\tilde{C}_{1}=u_{1}\rightarrow\cdots\rightarrow
u_{j}\rightarrow\cdots\rightarrow u_{p}\rightarrow\cdots\rightarrow u_{1}$ and
$\tilde{C}_{2}=u^{\prime}_{t}\rightarrow\cdots\rightarrow
u_{r}\rightarrow\cdots\rightarrow u_{q}\rightarrow u^{\prime}_{t}.$ Note that
$u_{1}\text{I}u_{j}\text{D}u^{\prime}_{t}$ (from the IDI path $\mathcal{P}),$
and $u_{1}$ is the initial point of $\tilde{C}$ (with respect to the standard
representation). So, $u_{1}<u^{\prime}_{t}<u_{j}.$ Again (after applying $f$
on $\tilde{\gamma})$ we have seen that $u_{j},u^{\prime}_{t}$ are the vertices
of two distinct cycles $\tilde{C}_{1},\tilde{C}_{2}$ respectively. Therefore,
$u^{\prime}_{t}$ is enclosed by the cycle $\tilde{C}_{1}.$ So,
$u^{\prime}_{t}$ is a singular point of the cycle $\tilde{C}_{1}.$ Also it is
clear that $u^{\prime}_{t}$ is the smallest singular point of the cycle
$\tilde{C}_{1}.$ As a result, $f(\tilde{\gamma})\in\mathcal{C}.$ In this case,
the number of cycles in $f(\tilde{\gamma})$ is increased by $1$ and
$w(\tilde{\gamma})=w(f(\tilde{\gamma})).$ now, it is easy to see that $f$ is a
weight preserving involution. Moreover, it is sign reversing as the number of
cycles is either increased or decreased by $1.$ This completes the proof. ∎

######  Example 8.

Here we explain the algorithm with one example. In Figure 2, the first diagram
has three cycles, $C_{1}=1\rightarrow 7\rightarrow 6\rightarrow 9\rightarrow
10\rightarrow 3\rightarrow 2\rightarrow 1,C_{2}=4\rightarrow 5\rightarrow
4\text{ and }C_{3}=8\rightarrow 11\rightarrow 8.$ As we observed above, in
$C_{1}$ the point of singularity are $6$ (as 1I7D6I9 is a IDI); $4,5$ (as
$4,5$ are vertices of the cycle $4\rightarrow 5\rightarrow 4\text{ and
}4,5,\text{ are enclosed by }C_{1});$ and $8$ (as $8$ is a vertex of the cycle
$8\rightarrow 11\rightarrow 8$ and $8$ is enclosed by the cycle $C_{1}).$
Also, the point of singularity of $C_{3}$ are $9,10$ (as $9,10$ are enclosed
by the cycle $8\rightarrow 11\rightarrow 8).$ But $C_{2}$ has no point of
singularity. Clearly, in the first diagram of Figure 2, $C_{1}$ is the
singular cycle with the smallest initial point (namely, $1$) and the smallest
point of singularity of $C_{1}$ is $4.$ So, according to the rule of $f,$ in
this example $C_{1},C_{2}$ plays the role of the cycle $C,C^{\prime}$
respectively, $v^{\prime}_{t}=4,v_{j}=10,v^{\prime}_{r}=5,v_{i}=3$ and the
cycle $C^{\prime\prime}=1\rightarrow 7\rightarrow 6\rightarrow 9\rightarrow
10\rightarrow 4\rightarrow 5\rightarrow 3\rightarrow 2\rightarrow 1.$

The second diagram in Figure 2 has two cycles namely,
$C^{\prime}_{1}=1\rightarrow 7\rightarrow 6\rightarrow 9\rightarrow
10\rightarrow 4\rightarrow 5\rightarrow 3\rightarrow 2\rightarrow 1\text{ and
}C^{\prime}_{2}=8\rightarrow 1\rightarrow 8.$ In $C^{\prime}_{1},$ the point
of singularity are $4$ (as $6\text{I}10\text{D}4\text{I}5$ is a IDI path); $6$
(as $1\text{I}7\text{D}6\text{I}10$ is a IDI path); $8$ (as $8$ is a vertex of
the cycle $8\rightarrow 10\rightarrow 8\text{ and }8\text{ is enclosed by
}C^{\prime}_{1}).$ The singular points of the cycle $C^{\prime}_{2}$ are $9$
and $10$ (as $9,10$ are the vertices of the cycle $C^{\prime}_{1}$ and $9,10$
are enclosed by the cycle $C^{\prime}_{2}).$ Clearly, in the second diagram of
Figure 2, $C^{\prime}_{1}$ is the singular cycle with the smallest initial
point (namely, $1$) and the smallest point of singularity of $C^{\prime}_{1}$
is $4.$ So, according to the rule of $f,$ here
$\tilde{C}=C^{\prime}_{1}=1\rightarrow 7\rightarrow 6\rightarrow 9\rightarrow
10\rightarrow 4\rightarrow 5\rightarrow 3\rightarrow 2\rightarrow 1$ and
$u^{\prime}_{t}=4.$ Therefore, if we apply $f$ on the first diagram of Figure
2, we get the second diagram of Figure 2 and conversely.

$\bf{1}$$\bf{2}$$\bf{3}$$\bf{4}$$\bf{5}$$\bf{6}$$\bf{7}$$\bf{8}$$\bf{9}$$\bf{10}$$\bf{11}$$\huge{\Updownarrow}$$\bf{1}$$\bf{2}$$\bf{3}$$\bf{4}$$\bf{5}$$\bf{6}$$\bf{7}$$\bf{8}$$\bf{9}$$\bf{10}$$\bf{11}$
Figure 2. After applying $f$ on the first diagram we get the second diagram
and conversely

### Acknowledgment

I would like to thank my mentor Prof. Arvind Ayyer for his constant support
and encouragement. The theorem is checked by “Sage”. We thank the authors for
generously releasing sage as an open-source package. The author was supported
by NBHM Post Doctoral Fellowship grant 0204/52/2019/RD-II/339. The author
sincerely thank the anonymous referees for their valuable suggestions and
comments.

## References

  * [1] A. Ayyer, _Determinants and perfect matchings_ , J. Combin. Theory Ser. A 120 (2013), no. 1, 304–314.
  * [2] S. Bera, _Enumeration of weighted paths on a digraph and block hook determinant_ , Spec. Matrices 9 (2021), 197–211.
  * [3] S. Bera and S. K. Mukherjee, _Combinatorial proofs of some determinantal identities_ , Linear Multilinear Algebra 66 (2018), no. 8, 1659–1667.
  * [4] G. Bhatnagar and C. Krattenthaler, _Spiral determinants_ , Linear Algebra Appl. 529 (2017), 374–390.
  * [5] A. R. Brualdi and D. Cvetkovic, _A combinatorial approach to matrix theory and its application_ , Discrete Mathematics and Its Applications, vol. 52, CRC Press, Boca Raton, London, New York, 2009.
  * [6] D. M. Jackson and I.P. Goulden, _A formal calculus for the enumerative system of sequences-I. combinatorial theorems_ , Studies in Applied Mathematics 61 (1979), 141–178.
  * [7] C. Krattenthaler, _Advanced determinant calculus_ , Sém. Lothar. Combin. 42 (1999), Art. B42q, 67.
  * [8] by same author, _Advanced determinant calculus: a complement_ , Linear Algebra Appl. 411 (2005), 68–166.
  * [9] P. D. Powell, _Calculating determinants of block matrices_ , arXiV:1112.4379v1 [math.RA] (2011).
  * [10] R. P. Stanley, _A matrix for counting paths in acyclic digraphs_ , J. Combin. Theory Ser. A 74 (1996), no. 1, 169–172.
  * [11] D. Zeilberger, _A combinatorial approach to matrix algebra_ , Discrete Mathematics 56 (1985), 61–72.
  * [12] D. Zheng, _Matrix methods for determinants of Pascal-like matrices_ , Linear Algebra Appl. 577 (2019), 94–113.