# Ultraviolet Freeze-in with a Time-dependent Inflaton Decay

Basabendu Barman  Nicolás Bernal
Yong Xu  and Óscar Zapata

###### Abstract

It is typically assumed that during reheating the inflaton decays with a
constant decay width. However, this is not guaranteed and can have a strong
impact on the dark matter (DM) genesis. In the context of the ultraviolet (UV)
freeze-in mechanism, if the operators connecting the dark and visible sectors
are of sufficiently high mass dimension, the bulk of the DM abundance is
produced during and not after reheating. We study here the impact of a time-
dependent decay width of the inflaton on the DM abundance, emphasizing the
differences with respect to the cases where the decay is either instantaneous
or constant. We also provide concrete examples for DM production via UV
freeze-in, e.g., from 2-to-2 scatterings of standard model particles, or from
inflaton scatterings or decays, elucidating how the time-dependence influences
the DM yield.

“All war is a symptom of man’s failure as a thinking animal.”
-John Steinbeck, ‘Once There Was A War’ (1958)

PI/UAN-2022-710FT

##  1 Introduction

The known baryonic matter can only explain $\sim 15\%$ of the total matter
budget of the Universe [1]; the rest is referred as dark matter (DM) [2, 3],
and is currently one of the most profound mysteries in both particle physics
and cosmology. Regarding the nature of DM, the weakly interacting massive
particles (WIMPs) are one of the most prominent candidates [4] (for a review,
see e.g. Refs. [5, 6, 7]). In the WIMP scenario, DM particles carry an
interaction strength at the electroweak scale, which allow them to thermalize
with the baryon-photon plasma in the early Universe and then freeze out,
reproducing the observed relic density $\Omega_{\rm DM}h^{2}\simeq 0.12$ [1].
This scenario is very appealing since several extensions of the standard model
(SM) of particle physics naturally include WIMPs. However, currently the
strong observational constraints on the typical WIMP parameter space motivate
quests beyond this paradigm (see, e.g., Refs. [8, 9]).

Alternatively to WIMPs, feebly interacting massive particles (FIMPs), which
couple to the SM sector very feebly evading therefore the current experimental
constraints, have been attracting heated discussion recently [10, 11, 12, 13,
14, 15]. In the early Universe, FIMPs can be generated from either the decay
or annihilation of states in the visible sector. When the SM temperature
becomes smaller than the typical mass scale of the interaction (i.e. the
maximum of the DM and the mediator mass), the generation process becomes
Boltzmann suppressed, giving rise to a constant comoving DM number density;
such a scenario is referred as the freeze-in mechanism [14].

The FIMP paradigm requires very suppressed interaction rates between the dark
and visible sectors, which can be achieved either in its infrared version via
small couplings (typically of the order $\sim 10^{-11}$), or in its
ultraviolet (UV) version via non-renormalizable operators [16], suppressed by
a high mass scale. The latter scenario is particularly interesting, as the DM
yield is sensitive to the highest temperature $T_{\text{max}}$ reached by the
SM plasma [17], controlled by the dynamics of the inflaton decay. In the
sudden decay approximation for the inflaton, $T_{\text{max}}$ corresponds to
the reheating temperature $T_{\text{rh}}$, characterizing the onset of the
radiation-dominance era. However, once away from the instantaneous reheating,
$T_{\text{max}}$ can be significantly larger than $T_{\text{rh}}$.

It is interesting to note that both $T_{\text{max}}$ and $T_{\text{rh}}$ are
controlled by the inflaton dynamics during reheating, in particular by its
dissipation rate $\Gamma_{\phi}$, typically assumed to be constant. However,
this should not be the case. In scenarios where the inflaton decays via
higher-order operators or oscillates around a potential steeper than
quadratic, its decay width may not be constant, but features a time dependence
[18, 19, 20, 21, 22, 23]. In particular, it was recently shown that in
contrast to the conventional case where during reheating the SM temperature
scales as $T(a)\propto a^{-3/8}$ (where $a$ corresponds to the scale factor),
the SM temperature could show a non-trivial dependence with the scale factor
[20].

Without loss of generality, we parametrize the dissipation rate as
$\Gamma_{\phi}(a,\,T)\propto a^{k}\,T^{q}$. Note that the conventional result
with constant $\Gamma_{\phi}$ corresponds to the case where $q=k=0$. In
general, one may expect a varying decay rate. For example, the dynamics of a
coherently oscillating scalar field in the early Universe can be affected by
the thermal environment due to thermal modification to the effective
potential, non-perturbative particle production or non-topological effects
[19]. If $\phi$ oscillates with a zero-temperature mass effective potential
$V\propto m_{\phi}^{2}\,\phi^{2}$, it is possible to have
$\Gamma_{\phi}\propto T$, in the limit the effective (temperature corrected)
mass $m_{\phi}^{\text{eff}}\sim m_{\phi}\ll T$ [19]. In these cases, it is
possible to realize $k=0$ and $q=1$. On the other hand, if the inflaton
oscillates in the vicinity of a potential steeper than quadratic, the decay
rate also typically features a scale-factor dependence due to the field
dependence of inflaton mass [21]. This can lead to a variety of values of $k$
depending on the spin of the decay products111Note that $\Gamma_{\phi}\propto
m_{\phi}$ or $1/m_{\phi}$ for fermions or bosons in the final state,
respectively. and the shape of the inflaton potential during reheating. For
example, for reheating in a quartic potential, one has $\Gamma_{\phi}\propto
a^{3/4}$ $(a^{-3/4})$ if the inflaton decays to bosons (fermions) [21]; this
corresponds to $q=0$ and $k=3/4$, or $q=0$ and $k=-3/4$. We consider a general
parametrization for the decay rate which captures a variety of dynamics during
reheating, including the effect of $i)$ decays via higher-order operators [20,
24], $ii)$ shapes of inflaton potentials during reheating [21, 25], as well as
$iii)$ feedback from the thermal background [19, 26].

The nontrivial behavior of the SM temperature during reheating leads to a
diverse DM phenomenology, in particular for those scenarios where its yield is
sensitive to the (highest) cosmic temperature after inflation. In particular,
if the operators connecting the dark and visible sectors are of sufficiently
high mass dimension, the bulk of the DM could be produced during and not after
reheating. This has been characterized by defining a “boost factor” $B$ for
the DM relic density, which is the ratio of the DM abundance taking into
account non-instantaneous reheating relative to the abundance in the sudden
decay approximation [27].222Another boost factor could be defined due to the
effects of thermalization and number-changing processes in the dark sector.
They can also have a strong impact, in particular enhancing the DM relic
abundance by several orders of magnitude [28]. For the sake of simplicity we
are assuming here no sizable self-interactions within the dark sector. It was
initially pointed out that such boosts only depend on the mass dimension of
the operator and the ratio $T_{\text{max}}/T_{\text{rh}}$, however the
equation-of-state parameter of the inflaton during reheating also plays a
major role [29, 30, 31, 32]. Subsequent papers have explored the impact of
this boost factor in specific models [33, 34, 35, 36, 37, 38, 39, 40, 41, 42,
43, 44, 45, 46, 47].

In this work we demonstrate that the phenomenology of the UV freeze-in
paradigm strongly depends on the dynamics of the inflaton decay during
reheating. Using a general parametrization for the inflaton decay width, we
carefully study the DM production during reheating and compute the
corresponding boost factors, by comparing it to the approximation where the
inflaton decays suddenly. In a completely model-independent fashion, we show
that the DM yield can have a power-law boost in the ratio
$T_{\text{max}}/T_{\text{rh}}$, particularly important if the SM temperature
drops fast enough during reheating. We typically concentrate on the situation
where the initial radiation density is negligibly small, and the total energy
density of the Universe is dominated by the inflaton during reheating. In such
a framework, we find that the nontrivial decay dynamics of the inflaton, and
consequently the boost factor, is completely determined by some combination of
the exponents $q$ and $k$, together with the two Hubble scales $H_{I}$ and
$H_{\text{rh}}$ defined at the beginning and the end of the reheating epoch,
respectively. Finally, we provide some realistic examples of DM production in
the early Universe, namely, gravitational UV freeze-in, inflaton scattering,
and inflaton decay, to demonstrate how a time-dependent inflaton decay width
can influence the DM yield in each of these cases.

The paper is organized as follows. In Sect. 2, we briefly revisit the
phenomenology of UV freeze-in where DM is produced out of annihilations of SM
particles, in the conventional approximation in which the inflaton suddenly
decays. In Sect. 3, we solve the Boltzmann equation to compute the background
evolution in the presence of a time-dependent inflaton decay rate. The
phenomenology of the modified FIMP DM yield is studied in Sect. 4. In Sec. 5
we provide three concrete instances to exemplify our results. Finally, we
summarize our findings in Sect. 6.

##  2 UV Freeze-in in the Sudden Decay Approximation

The evolution of the DM number density $n_{\text{DM}}$ is governed by a
Boltzmann equation (BEQ), that can be written in a generalized form as

$\frac{dn_{\text{DM}}}{dt}+3\,H\,n_{\text{DM}}=\gamma(T)\,,$ (2.1)

where $\gamma(T)$ corresponds to the DM production rate density out of SM
particles, as a function of the bath temperature $T$. In case where the
Universe energy density is dominated by SM radiation, the Hubble expansion
rate $H$ is given by

$H(T)=\sqrt{\frac{\rho_{R}(T)}{3\,M_{P}^{2}}}\,,$ (2.2)

with the SM radiation energy density

$\rho_{R}(T)=\frac{\pi^{2}}{30}\,g_{\star}(T)\,T^{4}\,,$ (2.3)

where $g_{\star}$ is the number of relativistic degrees of freedom
contributing to the SM energy density [48], and $M_{P}$ the reduced Planck
mass. In the case where the SM entropy is conserved, it is instructive to
solve Eq. (2.1) in terms of the DM yield $Y\equiv n_{\text{DM}}/s$, defined as
a ratio of DM number to entropy density of the Universe, where
$s(T)\equiv\frac{2\pi^{2}}{45}g_{\star s}\,T^{3}$, with $g_{\star s}(T)$ being
the number of relativistic degrees of freedom contributing to the SM entropy
[48]. At temperatures much higher than the electroweak scale,
$g_{\star}=g_{\star s}=106.75$. On substituting in Eq. (2.1) one obtains

$\frac{dY}{dT}=-\frac{\gamma(T)}{H(T)\,T\,s(T)}\,.$ (2.4)

Now, in case of UV freeze-in, the temperature $T$ of the thermal bath is
always lower than the typical mass scale of the interaction, implying that the
operators through which the DM communicates with the visible sector are non-
renormalizable. In such a scenario, the reaction rate for DM produced from the
SM bath can be parametrized as [16]

$\gamma(T)=\frac{T^{n}}{\Lambda^{n-4}}\,,$ (2.5)

where $\Lambda$ is a dimensionful parameter which is the effective interaction
scale between the DM and the SM, coming from an operator of mass dimension
$d=2+n/2$, (with $n\geq 6$).333Detection prospects of UV freeze-in via non-
renormalizable operators in an early matter dominated era have been discussed
in Ref. [49]. It is thus implied that the effective description is valid at
temperatures $T<\Lambda$.

In the approximation where the inflaton decays instantaneously, the reheating
temperature $T_{\text{rh}}$ corresponds to the maximum temperature achieved by
the thermal bath and characterizes the onset of the radiation-dominated era.
For SM temperatures much smaller than $T_{\text{rh}}$ (i.e., $T\ll
T_{\text{rh}}$), it is possible to solve Eq. (2.4) analytically, obtaining

$Y_{0}\equiv Y(T\ll T_{\text{rh}})=\frac{135}{2\pi^{3}(n-5)\,g_{\star
s}}\sqrt{\frac{10}{g_{\star}}}\,\frac{M_{P}\,T_{\text{rh}}^{n-5}}{\Lambda^{n-4}}\,,$
(2.6)

where $Y_{0}$ denotes the DM yield at the end of reheating and assuming
$m_{\text{DM}}\ll T_{\text{rh}}$, $m_{\text{DM}}$ being the DM mass. It is
important to emphasize that we have also assumed a DM abundance initially
negligible and a production cross section sufficiently small such that DM
remains out of chemical equilibrium with the SM bath. Finally, to match the
whole observed abundance, the DM yield has to be fixed so that
$m_{\text{DM}}\,Y_{0}=\Omega_{\text{DM}}h^{2}\frac{1}{s_{0}}\,\frac{\rho_{c}}{h^{2}}\simeq
4.3\times 10^{-10}$ GeV, where $m_{\text{DM}}$ is the DM mass in GeV,
$\rho_{c}\simeq 1.1\times 10^{-5}h^{2}$ GeV/cm3 is the critical energy
density, $s_{0}\simeq 2.9\times 10^{3}$ cm-3 is the entropy density at
present, and $\Omega_{\text{DM}}h^{2}\simeq 0.12$ [1].

##  3 Non-instantaneous Reheating with a Time-dependent Inflaton Decay

Although reheating is commonly approximated to be sudden, in reality, the
decay of the inflaton is non-instantaneous and typically characterized by an
exponential decay law [50]. Therefore, the evolution of inflaton and radiation
energy densities ($\rho_{\phi}$ and $\rho_{R}$, respectively) is governed by a
set of coupled BEQs that read [17]444We would like to mention that we have
focused on a scenario where inflaton decays perturbatively. Potentially
relevant non-perturbative preheating processes [51] will not be studied here.

$\displaystyle\frac{d\rho_{\phi}}{dt}+3\,H\,\rho_{\phi}=-\Gamma_{\phi}\,\rho_{\phi}\,,$
(3.1)
$\displaystyle\frac{d\rho_{R}}{dt}+4\,H\,\rho_{R}=+\Gamma_{\phi}\,\rho_{\phi}\,,$
(3.2)

with, in this case, the Hubble expansion rate $H$ is given by

$H^{2}=\frac{\rho_{R}+\rho_{\phi}}{3\,M_{P}^{2}}\,,$ (3.3)

and $\Gamma_{\phi}$ being the total perturbative decay width of the inflaton.
It is also assumed that during reheating (i.e., for $a_{I}<a<a_{\text{rh}}$,
where $a_{I}$ and $a_{\text{rh}}$ correspond to the scale factors at the end
of inflation and at the end of reheating, respectively), the inflaton energy
density scales as non-relativistic matter, i.e., $\rho_{\phi}(a)\propto
a^{-3}$. This scaling is characteristic of an inflaton oscillating in a
quadratic potential. Therefore, during reheating

$\rho_{\phi}(a)\simeq 3M_{P}^{2}\,H_{I}^{2}\left(\frac{a_{I}}{a}\right)^{3},$
(3.4)

with $H_{I}\equiv H(a=a_{I})$ denoting the Hubble parameter at the end of
inflation. As advocated earlier, generalizing the widely used assumption of
constant $\Gamma_{\phi}$, here we parametrize the inflaton decay as a function
of temperature and scale factor as

$\Gamma_{\phi}(a,\,T)=\mathcal{C}\left(\frac{a}{a_{\text{rh}}}\right)^{k}\left(\frac{T}{T_{\text{rh}}}\right)^{q}H_{\text{rh}}\,,$
(3.5)

with $\mathcal{C}=\mathcal{C}(k,\,q)$ a constant that will be conveniently
fixed in the following.555Note that this parametrization slightly differs from
the one in Ref. [20]. We would like to stress that the previous equation is
just a convenient parametrization, the exact expression could be computed once
a full model is fixed. For instance, the present parametrization is not valid
before the thermalization of the SM bath. The decay width can be recasted in
terms of the radiation energy density following Eq. (2.3) as

$\Gamma_{\phi}(a)=\mathcal{C}\left(\frac{30}{\pi^{2}\,g_{\star}}\right)^{q/4}\left(\frac{a}{a_{\text{rh}}}\right)^{k}\left(\frac{\rho_{R}^{1/4}(a)}{T_{\text{rh}}}\right)^{q}H_{\text{rh}}\,,$
(3.6)

where $q$ and $k$ are parameters given by the specific coupling of the
inflaton to SM particles. Note that if $k=q=0$, one actually reproduces the
instantaneous reheating scenario with
$\Gamma_{\phi}(T_{\text{rh}})=\mathcal{C}\,H_{\text{rh}}$. In that sense, the
parameterization of $\Gamma_{\phi}$ in Eq. (3.5) is generic and captures both
the standard and some potentially unusual features of the reheating dynamics.
The reheating period ends at the so-called reheat temperature $T_{\text{rh}}$,
which can be defined using $\rho_{R}(T_{\text{rh}})\equiv
3M_{P}^{2}\,H_{\text{rh}}^{2}$,

$T_{\text{rh}}^{2}=\frac{3}{\pi}\sqrt{\frac{10}{g_{\star}}}\,M_{P}\,H_{\text{rh}}\,,$
(3.7)

which also marks the onset of the radiation-domination era.

During reheating, it is possible to analytically solve Eqs. (3.1) and (3.2) in
order to extract the radiation energy density as

$\displaystyle\rho_{R}(a)$
$\displaystyle=\Biggl{[}\left(3\,M_{P}^{2}\,H_{\text{rh}}^{2}\right)^{\frac{4-q}{4}}\frac{\mathcal{C}}{2}\frac{4-q}{5+2k-2q}\,\frac{H_{I}}{H_{\text{rh}}}\left(\frac{a_{I}}{a_{\text{rh}}}\right)^{k}\left(\Bigl{(}\frac{a}{a_{I}}\Bigr{)}^{\frac{5+2k-2q}{2}}-1\right)$
$\displaystyle\qquad+\rho_{R}(a_{I})^{\frac{4-q}{4}}\Biggr{]}^{\frac{4}{4-q}}\,\left(\frac{a_{I}}{a}\right)^{4},$
(3.8)

where it was assumed that the total energy density is dominated by inflatons
during reheating. As just after inflation one expects a completely subdominant
SM radiation energy density, that is $\rho_{R}(a_{I})\simeq 0$, here we will
focus on the case $q<4$. For completeness, however, we would like to mention
the consequences of having $q>4$. As discussed in Ref. [20], although the
$q>4$ scenario is governed by the same Eq. (3), it gives rise to a very
different outcome. This is typically due to the fact that unlike the $q<4$
scenario, the initial radiation energy density is no more negligible here.
Now, for $5-2q+2k>0$, the temperature continues to decrease during reheating,
with a sharp increase towards the end of reheating showing a sudden completion
of dissipation. On the other hand, with $5-2q+2k<0$, the dissipation is never
complete and one ends up with the usual instantaneous reheating. Now, on top
of $q<4$, we additionally require $5-2q+2k>0$ to guarantee an efficient energy
transfer from the inflaton to the SM radiation. Under such considerations one
has

$\rho_{R}(a)\simeq
3\,M_{P}^{2}\,H_{\text{rh}}^{2}\left[\frac{\mathcal{C}}{2}\frac{4-q}{5-2q+2k}\,\frac{H_{I}}{H_{\text{rh}}}\left(\frac{a_{I}}{a_{\text{rh}}}\right)^{k}\right]^{\frac{4}{4-q}}\left(\frac{a_{I}}{a}\right)^{\frac{6-4k}{4-q}}.$
(3.9)

At the end of reheating (i.e., at $a=a_{\text{rh}}$), the inflaton and the
radiation energy densities are equal,
$\rho_{R}(a_{\text{rh}})=\rho_{\phi}(a_{\text{rh}})=3M_{P}^{2}\,H_{\text{rh}}^{2}$,
and therefore Eq. (3.9) implies

$\frac{a_{I}}{a_{\text{rh}}}=\left(\frac{2}{\mathcal{C}}\,\frac{5+2k-2q}{(4-q)}\,\frac{H_{\text{rh}}}{H_{I}}\right)^{\frac{2}{3}}\equiv\left(\frac{H_{\text{rh}}}{H_{I}}\right)^{\frac{2}{3}},$
(3.10)

where in the last term we have used the Hubble-scale factor relation during
reheating. The constant $\mathcal{C}$ is thus fixed from now on as

$\mathcal{C}=\frac{2\,\left(5-2q+2k\right)}{4-q}\,,$ (3.11)

which in turn allows to rewrite Eq. (3.9) as

$\rho_{R}(a)\simeq
3\,M_{P}^{2}\,H_{\text{rh}}^{2}\left(\frac{a_{\text{rh}}}{a}\right)^{\frac{2\,(3-2k)}{4-q}}.$
(3.12)

It is worth mentioning that for $k>3/2$, $\rho_{R}$ and therefore the SM
temperature continuously increases during reheating, and therefore
$T_{\text{rh}}$ corresponds to the maximum temperature reached by the SM
plasma. Additionally, if $k=3/2$, the SM temperature remains constant during
the inflaton-dominated era. Finally, note that for $k<3/2$, in the first
stages of reheating, the SM bath reaches a maximal temperature of
$T_{\text{max}}$ [17], given by

$T_{\text{max}}=T_{\text{rh}}\left(\frac{H_{I}}{H_{\text{rh}}}\right)^{\frac{3-2k}{3\,(4-q)}},$
(3.13)

which is higher than $T_{\text{rh}}$.

It is interesting to note that the decay width could be simplified to the form

$\Gamma_{\phi}(a)\simeq\frac{5-2x}{2}\,H_{\text{rh}}\left(\frac{a_{\text{rh}}}{a}\right)^{x},$
(3.14)

where

$x\equiv\frac{3q-8k}{2\,(4-q)}\,,$ (3.15)

with $x<5/2$ coming from the condition $5-2q+2k>0$. Additionally, one still
has to guarantee that $q<4$. We emphasize that a constant decay width is
reached in the case $k=q=0$, but also more generally for $x=0$, corresponding
to $k=\frac{3}{8}\,q$. Cases $x<0$ and $x>0$ generate an inflaton decay width
that increases or decreases with time, respectively. The SM radiation energy
density becomes, in turn,

$\rho_{R}(a)\simeq
3\,M_{P}^{2}\,H_{\text{rh}}^{2}\left(\frac{a_{\text{rh}}}{a}\right)^{\frac{3+2x}{2}}.$
(3.16)

Before moving on, we would like to emphasize that, as expected, the present
scenario demands a set of three free parameters that could conveniently be
chosen to be $x$, together with the two Hubble scales: $H_{I}$ and
$H_{\text{rh}}$. The Planck collaboration puts an upper limit on the
inflationary scale $H_{I}\leq 2.5\times 10^{-5}~{}M_{P}$ [52] from the non-
observation of CMB tensor modes. On the other hand, by combining light element
abundance measurements with CMB and large-scale structure data, a fairly
robust lower limit on the reheating temperature $T_{\text{rh}}\gtrsim 4$ MeV
can be set at 95% CL [53, 54, 55] that in turn puts a lower bound on
$H_{\text{rh}}\gtrsim 3.0\times 10^{-42}~{}M_{P}$. In the upper panels of Fig.
1 the evolution of the inflaton (blue) and SM radiation (black) energy
densities as a function of the scale factor, for $x=-3$ (left), $x=-3/2$
(center) and $x=0$ (right). We take $H_{I}=10^{8}$ GeV and
$H_{\text{rh}}=10^{3}$ GeV. The lower panels show the corresponding evolution
of the photon temperature. The thick solid lines show the numerical solutions
of the system of Boltzmann equations, whereas the red dotted lines the
analytical approximations.

Figure 1: Upper panels: Evolution of the inflaton (blue) and SM radiation
(black) energy densities as a function of the scale factor, for $x=-3$ (left),
$x=-3/2$ (center) and $x=0$ (right), taking $H_{I}=10^{8}$ GeV and
$H_{\text{rh}}=10^{3}$ GeV. Lower panels: The corresponding evolution of the
photon temperature.

##  4 UV Freeze-in with a Time-dependent Inflaton Decay

In this section, we address the effect of a generic time-dependent inflaton
decay rate on the DM production via UV freeze-in.666We emphasize that here we
are focusing on a non-instantaneous decay of the inflaton, but on an
instantaneous thermalization of the SM plasma [56, 57, 58]. The SM particles
do not necessarily thermalize suddenly, and thus they could be initially
distributed with smaller occupation numbers and harder momenta [59, 60, 61,
62]. Since the SM entropy is not conserved during reheating, it is more
convenient to rewrite Eq. (2.1) in term of the comoving DM number density
$N\equiv n_{\text{DM}}\times a^{3}$ as777We note that DM could also be
produced via freeze-in from parametric resonances (i.e. preheating) [63].

$\frac{dN}{da}=\frac{a^{2}\gamma}{H}\simeq\left(\frac{a_{\text{rh}}^{3\,n/n_{c}}}{a_{I}}\right)^{3/2}\frac{T_{\text{rh}}^{n}}{H_{I}\,\Lambda^{n-4}}\,a^{\frac{7}{2}-\frac{9}{2}\frac{n}{n_{c}}},$
(4.1)

with

$n_{c}\equiv\frac{36}{3+2x}=9\,\frac{4-q}{3-2k}\,,$ (4.2)

defined for $x\neq-3/2$ (or equivalently, for $k\neq 3/2$). The case $x=-3/2$,
corresponding to a constant SM temperature during reheating, will be analyzed
later on.

Now, the DM number at the end of reheating can be analytically obtained by
integrating Eq. (4.1) over the range $a_{I}<a<a_{\text{rh}}$,

$N(a_{\text{rh}})\simeq\begin{dcases}\frac{2}{9}\frac{n_{c}}{n_{c}-n}\,\frac{a_{\text{rh}}^{3}\,T_{\text{rh}}^{n}}{H_{I}\,\Lambda^{n-4}}\left(\frac{a_{\text{rh}}}{a_{I}}\right)^{\frac{3}{2}}\left[1-\left(\frac{a_{I}}{a_{\text{rh}}}\right)^{\frac{9}{2}\frac{n_{c}-n}{n_{c}}}\right]&\text{
for }n\neq n_{c}\,,\\\
\frac{a_{\text{rh}}^{3}\,T_{\text{rh}}^{n_{c}}}{H_{I}\,\Lambda^{n_{c}-4}}\,\left(\frac{a_{\text{rh}}}{a_{I}}\right)^{\frac{3}{2}}\ln\frac{a_{\text{rh}}}{a_{I}}&\text{
for }n=n_{c}\,,\end{dcases}$ (4.3)

where a DM lighter than $T_{\text{rh}}$ was assumed. The case where DM is
heavier than $T_{\text{rh}}$, but still lighter than $T_{\text{max}}$, will be
considered separately. It is important to note that we are also taking a
subdominant initial DM number density, $N(a_{I})\simeq 0$, a common assumption
at the end of inflation. Although the SM entropy density is not conserved when
the inflaton is decaying, one can further define the DM yield at
$a=a_{\text{rh}}$, namely
$Y(a_{\text{rh}})=n_{\text{DM}}(a_{\text{rh}})/s(T_{\text{rh}})=N(a_{\text{rh}})/(a_{\text{rh}}^{3}\,s(T_{\text{rh}}))$
as

$Y(a_{\text{rh}})\simeq\begin{dcases}\frac{n_{c}}{n-n_{c}}\,\frac{5}{\pi^{2}\,g_{\star
s}}\,\frac{T_{\text{rh}}^{n-3}}{H_{\text{rh}}\,\Lambda^{n-4}}\,\Biggl{[}\left(\frac{H_{I}}{H_{\text{rh}}}\right)^{3\,\frac{n-n_{c}}{n_{c}}}-1\Biggr{]}&\text{
for }n\neq n_{c}\,,\\\ \frac{45}{(2\,\pi^{2})\,g_{\star
s}}\,\frac{T_{\text{rh}}^{n_{c}}}{T_{\text{rh}}^{3}\,H_{\text{rh}}\,\Lambda^{n_{c}-4}}\,\ln\left(\frac{H_{I}}{H_{\text{rh}}}\right)^{\frac{2}{3}}\,&\text{
for }n=n_{c}\,.\end{dcases}$ (4.4)

After the end of reheating, the SM entropy is conserved and therefore
$Y(a_{\text{rh}})$ remains constant.

To quantify the DM production during reheating, we define a boost factor $B$
for the DM relic density, which is the ratio of the DM abundance taking into
account non-instantaneous reheating relative to the abundance under the sudden
decay approximation [27, 29]. Therefore, the total DM yield $Y_{\text{total}}$
corresponds to the production after and during reheating, so that

$Y_{\text{total}}=Y_{0}+Y(a_{\text{rh}})=Y_{0}\times\left(1+B\right)\,.$ (4.5)

By comparing Eqs. (2.6) with (4.4), the boost factor is estimated to be

$B\equiv\frac{Y(a_{\text{rh}})}{Y_{0}}\simeq\begin{dcases}\frac{2}{9}\,\frac{(n-5)\,n_{c}}{n-n_{c}}\left(\frac{H_{I}}{H_{\text{rh}}}\right)^{3\,\frac{n-n_{c}}{n_{c}}}&\text{
for }n>n_{c}\,,\\\
\frac{2}{3}\,(n_{c}-5)\,\ln\left(\frac{H_{I}}{H_{\text{rh}}}\right)\,&\text{
for }n=n_{c}\,,\\\ \frac{2}{9}\,\frac{(n-5)\,n_{c}}{n_{c}-n}&\text{ for
}n<n_{c}\,,\end{dcases}$ (4.6)

which tells that a large boost (i.e., a power-law enhancement in the ratio
$H_{I}/H_{\text{rh}}$) in the DM production during reheating appears for the
case where $n>n_{c}$. If $n=n_{c}$, the boost is logarithmic in the ratio
$H_{I}/H_{\text{rh}}$, whereas for $n<n_{c}$ an $\mathcal{O}(1)$ boost is
expected. We note that such a boost depends on $x$, $n$, and the ratio of the
Hubble expansion rates, but not on $\Lambda$ or $m_{\text{DM}}$. The case
$x=0$ (i.e., an inflaton with an energy density that during reheating scales
like non-relativistic matter and decays with a constant decay width) matches
with previous results reported in the literature [27, 29]. Note that in the
case with $x>-3/2$ (or equivalently $k<3/2$), the SM energy density
monotonically decreases during reheating, so there is a well defined
$T_{\text{max}}$, given by Eq. (3.13). In this case, one can further write the
boost factor in Eq. (4.6) in terms of $T_{\text{max}}/T_{\text{rh}}$ as

$B\simeq\begin{dcases}\frac{2}{9}\,\frac{(n-5)\,n_{c}}{n-n_{c}}\left(\frac{T_{\text{max}}}{T_{\text{rh}}}\right)^{n-n_{c}}&\text{
for }n>n_{c}\,,\\\
\frac{2}{9}\,n_{c}\,(n_{c}-5)\,\ln\left(\frac{T_{\text{max}}}{T_{\text{rh}}}\right)\,&\text{
for }n=n_{c}\,,\\\ \frac{2}{9}\,\frac{(n-5)\,n_{c}}{n_{c}-n}&\text{ for
}n<n_{c}\,,\end{dcases}$ (4.7)

where
$H_{I}/H_{\text{rh}}=\left(T_{\text{max}}/T_{\text{rh}}\right)^{n_{c}/3}$ has
been utilized.

We now turn to the special case $x=-3/2$ (equivalently to $k=3/2$). As evident
from Eq. (3.16), this corresponds to a constant SM temperature during
reheating. Thus, $x=-3/2$ marks the transition between the cases where the SM
energy density grows with time during reheating $(x<-3/2)$, and where there is
no period of increasing temperature $(x>-3/2)$. In this case one can write the
evolution of DM number as

$\frac{dN}{da}\simeq
a_{I}^{-3/2}\frac{T_{\text{rh}}^{n}}{H_{I}\,\Lambda^{n-4}}\,a^{\frac{7}{2}},$
(4.8)

following Eq. (4.1), which eventually leads to a DM number

$N(a_{\text{rh}})\simeq\frac{2}{9}\frac{T_{\text{rh}}^{n}}{H_{I}\,\Lambda^{n-4}}\,a_{I}^{3}\left[\left(\frac{a_{\text{rh}}}{a_{I}}\right)^{\frac{9}{2}}-1\right],$
(4.9)

at the end of reheating, again assuming $m_{\text{DM}}<T_{\text{rh}}$.
Following the same prescription as before, the DM yield can be obtained as

$Y(a_{\text{rh}})\simeq\frac{5}{\pi^{2}\,g_{\star
s}}\,\frac{T_{\text{rh}}^{n-3}}{H_{\text{rh}}\,\Lambda^{n-4}}\,.$ (4.10)

Finally, the boost factor in this case reads

$B\simeq\frac{2}{9}\,(n-5)\,,$ (4.11)

and is always $\mathcal{O}(1)$. This can be understood by noticing that, even
if DM is constantly produced during reheating and at the same rate (the SM
temperature being a constant), the DM yield suffers from the entropy dilution
from the decay of the inflaton. Therefore, the production is largely dominated
by late times, close to the end of reheating, and hence a small boost is
expected.

In Fig. 2 we depict the parameter space of our interest in the $[k,\,q]$
plane, where the black-slanted straight contours correspond to
$n_{c}=\\{8,\,10,\,12,\,14\\}$ from bottom to top. The shaded gray region is
excluded from the requirement of $q-4<0$ and $5-2q+2k>0$, which corresponds to
a vanishing radiation energy density just after inflation. Let us remind that
the parameter space below (above) the line $n_{c}=12$, or equivalently $x=0$,
is associated with an inflaton decay width that increases (decreases) with
time. Since, as mentioned before, a power-law enhanced boost in the DM yield
is achievable only in the case where $n_{c}<n$, hence below these contours a
sizable boost factor is viable for operators with a given $n>n_{c}$. Finally,
we note that, as expected, a large boost can only appear in the case where
$k<3/2$, where the SM bath could reach a temperature $T_{\text{max}}$ much
larger than $T_{\text{rh}}$.

Figure 2: In the $[k,\,q]$ plane the black straight contours (cf. Eq. (4.2))
correspond to different choices of $n_{c}:\\{8\,,10\,,12\,,14\\}$ from bottom
right to top left. The shaded regions are disallowed from the requirement of
having $q<4$ and $5-2q+2k>0$. Below the red dotted straight horizontal line at
$k=3/2$ temperature of the Universe decreases during reheating, while above,
the temperature increases throughout reheating (see text).

Before closing this section, we focus on the case where $T_{\text{max}}\gg
m_{\text{DM}}\gg T_{\text{rh}}$, which can happen only for $x>-3/2$ (or
equivalently $k<3/2$). In this scenario, DM is produced during, but not after
reheating. The DM number could be estimated by integrating Eq. (4.1) from
$a_{I}$ to $a_{\text{DM}}$, where $a_{\text{DM}}\equiv a(T=m_{\text{DM}})$ is
given by

$a_{\text{DM}}=a_{\text{rh}}\left(\frac{T_{\text{rh}}}{m_{\text{DM}}}\right)^{\frac{2n_{c}}{9}}.$
(4.12)

Thus, one gets

$N(a_{\text{DM}})\simeq\begin{dcases}\frac{2}{9}\frac{n_{c}}{n_{c}-n}\frac{a_{I}^{3}\,T_{\text{max}}^{n}}{H_{I}\,\Lambda^{n-4}}\left[\left(\frac{T_{\text{max}}}{m_{\text{DM}}}\right)^{n_{c}-n}-1\right]&\text{
for }n\neq n_{c}\,,\\\
\frac{2}{3}\frac{a_{\text{rh}}^{3}\,T_{\text{rh}}^{n_{c}}}{H_{\text{rh}}\,\Lambda^{n_{c}-4}}\,\ln\left[\left(\frac{T_{\text{rh}}}{m_{\text{DM}}}\right)^{\frac{n_{c}}{3}}\frac{H_{I}}{H_{\text{rh}}}\right]&\text{
for }n=n_{c}\,.\end{dcases}$ (4.13)

Hence the DM yield
$Y(a_{\text{DM}})=N(a_{\text{DM}})/(a_{\text{DM}}^{3}\,s(m_{\text{DM}}))$
reads

$Y(a_{\text{DM}})\simeq\begin{dcases}\frac{5}{\pi^{2}\,g_{\star
s}}\,\frac{n_{c}}{n_{c}-n}\,\frac{m_{\text{DM}}^{n-3}}{H_{I}\,\Lambda^{n-4}}\left(\frac{T_{\text{max}}}{m_{\text{DM}}}\right)^{\frac{n_{c}}{3}}&\text{
for }n<n_{c}\,,\\\ \frac{5}{\pi^{2}\,g_{\star
s}}\,\frac{n_{c}}{n-n_{c}}\frac{m_{\text{DM}}^{n-3}}{H_{\text{rh}}\,\Lambda^{n-4}}\left(\frac{T_{\text{max}}}{m_{\text{DM}}}\right)^{n-n_{c}}\left(\frac{T_{\text{rh}}}{m_{\text{DM}}}\right)^{\frac{n_{c}}{3}}&\text{
for }n>n_{c}\,.\end{dcases}$ (4.14)

It is important to note here that the DM yield is conserved after reheating
($a>a_{\text{rh}}$), but not during reheating (and in particular within the
range $a_{\text{DM}}<a<a_{\text{rh}}$), due to the entropy dilution from the
inflaton decay. Then, the DM yield at the end of reheating is

$Y(a_{\text{rh}})=Y(a_{\text{DM}})\,\frac{S(a_{\text{DM}})}{S(a_{\text{rh}})}=Y(a_{\text{DM}})\,\left(\frac{m_{\text{DM}}}{T_{\text{rh}}}\right)^{3}\left(\frac{T_{\text{rh}}}{m_{\text{DM}}}\right)^{\frac{2n_{c}}{3}}.$
(4.15)

We emphasize that in this case the DM is only produced during (and not after)
reheating, and therefore the boost factor can not be defined.

##  5 Scenarios with Time-dependent Inflaton Decay

Until now, we have discussed the consequences of time-dependent inflaton
decay, especially focusing on how it can modify the DM abundance completely
model-agnostic, using the generic parameterization for the inflaton decay
width (cf. Eq. (3.5)), together with a generic DM production rate in Eq.
(2.5). In this section, we would motivate how such DM production rates can
originate from a physically realizable framework. We focus on the following
three situations, where DM can be produced via $i)$ 2-to-2 scattering of the
bath particles through an $s$-channel graviton exchange, $ii)$ 2-to-2
scattering of the inflatons during reheating via a graviton exchange, and
$iii)$ decay of the inflaton field into DM final states. The first two cases
belong to pure gravitational production, where the DM particle production
occurs via gravitational interaction which is described in terms of coupling
of the energy-momentum tensor to the metric perturbation, where the latter is
identified as the quantum field for the spin-2 massless graviton [64, 65, 66].
This subsequently gives rise to interactions of all matter fields with the
graviton and also pair production of DM particles, satisfying the whole
observed DM abundance [67, 68, 69, 70, 71, 72, 73]. On the other hand, DM can
also be produced directly from inflaton decay with a tiny branching fraction,
if the inflaton couples to DM fields.888Even if the DM does not have a direct
coupling to the inflaton, one cannot ignore the radiative decay of inflaton
into DM [37]. In the following subsections, we shall see how a time-dependent
inflaton decay width influences DM production within the frameworks mentioned
above.

###  5.1 Gravitational UV Freeze-in

We examine the first case where the DM is gravitationally produced via the UV
freeze-in mechanism of 2-to-2 annihilation of SM particles mediated by the
$s$-channel exchange of gravitons. The interaction rate density for such a
process reads [67, 70, 69, 34, 74, 75, 73]

$\gamma=\alpha\,\frac{T^{8}}{M_{P}^{4}}\,,$ (5.1)

with $\alpha\simeq 1.9\times 10^{-4}$ (real scalar DM), $\alpha\simeq
1.1\times 10^{-3}$ (Dirac DM) or $\alpha\simeq 2.3\times 10^{-3}$ (vector DM).
Compared to Eq. (2.5), we find that this situation corresponds to $n=8$.999UV
freeze-in via SM scattering through the dilaton portal can give rise to $n=12$
for fermionic DM [76]. Thus, comparing with the general form of the boost
factor in Eq. (4.6), we see that the DM yield is boosted if $n_{c}<n=8$, which
corresponds to $x>3/4$. This is reflected in Fig. 3, where contours of
different choices of the boost factor are shown, for $n=8$. The left panel
describes the required ratio of Hubble parameter, namely $H_{I}/H_{\text{rh}}$
in order to yield a fixed boost factor $B$, while the right panel depicts the
same with respect to $T_{\text{max}}/T_{\text{rh}}$. Note that these two
ratios are connected by
$H_{I}/H_{\text{rh}}=\left(T_{\text{max}}/T_{\text{rh}}\right)^{n_{c}/3}$, cf.
Eq. (3.13). Small $\mathcal{O}(1)$ boost appear when $x<3/4$, whereas larger
boosts can happen in the opposite case. This figure was produced using the
analytical approximation in Eq. (4.6); we have checked the good agreement with
the fully numerical result.

Figure 3: Contours of boost factor $B$ for the gravitational UV freeze-in
scenario with SM scattering (i.e., $n=8$), where $B=10^{0}$, $10^{1}$,
$10^{2}$ and $10^{3}$ from bottom left to top right. Power-law boosts happen
for $3/4<x<5/2$. The gray-shaded area corresponding to $x>5/2$ cannot be
explored in the case of a negligible DM population at the end of inflation.
The left and right panel correspond to the ratios $H_{I}/H_{\text{rh}}$ and
$T_{\text{max}}/T_{\text{rh}}$, respectively.

###  5.2 Gravitational Inflaton Scattering

Similarly to SM scattering, during reheating, the entire observed DM abundance
can be generated by 2-to-2 annihilations of inflatons, mediated by the
$s$-channel exchange of gravitons. This corresponds to the second case in our
list of examples. The interaction rate density for DM production out of non-
relativistic inflatons then reads [71, 73, 72, 77, 78, 79]

$\gamma=\frac{\rho_{\phi}^{2}}{512\pi\,M_{P}^{4}}\,f\left(\frac{m_{\text{DM}}}{m_{\phi}}\right),$
(5.2)

where $m_{\phi}$ is the inflaton mass and

$f(y)\equiv\begin{dcases}\left(y^{2}+2\right)^{2}\sqrt{1-y^{2}}&\text{for real
scalars},\\\ y^{2}\left(1-y^{2}\right)^{3/2}&\text{for Dirac fermions},\\\
\frac{1}{8}\sqrt{1-y^{2}}\left(4+4y^{2}+19y^{4}\right)&\text{for vector
bosons.}\end{dcases}$ (5.3)

By writing Eq. (3.4) as a function $n_{c}$

$\rho_{\phi}(T)=3\,M_{P}^{2}\,H_{\text{rh}}^{2}\left(\frac{T}{T_{\text{rh}}}\right)^{\frac{2}{3}n_{c}},$
(5.4)

one can express the reaction density as

$\gamma=\frac{9}{512\,\pi}\,H_{\text{rh}}^{4}\left(\frac{T}{T_{\text{rh}}}\right)^{\frac{4}{3}n_{c}}f\left(\frac{m_{\text{DM}}}{m_{\phi}}\right)=\frac{9}{512\,\pi}\,H_{\text{rh}}^{4}\left(\frac{a_{\text{rh}}}{a}\right)^{6}f\left(\frac{m_{\text{DM}}}{m_{\phi}}\right).$
(5.5)

Again, comparing with Eq. (2.5), we see that in this case,
$\Lambda^{\frac{4\,(3-n_{c})}{3}}=\frac{9}{512\pi}\left(H_{\text{rh}}\,T_{\text{rh}}^{-\frac{n_{c}}{3}}\right)^{4}f$
and $n=\frac{4}{3}n_{c}$. Following the same methodology as before, we can
determine the DM yield at the end of reheating by solving the BEQ (namely Eq.
(4.1)), that reads

$Y(a_{\text{rh}})=\frac{135}{512\,\pi^{3}\,g_{\star
s}}\,\frac{H_{I}\,H_{\text{rh}}^{2}}{T_{\text{rh}}^{3}}\,f\left(\frac{m_{\text{DM}}}{m_{\phi}}\right),$
(5.6)

which is, interestingly, independent of $x$. This is expected from Eq. (5.5),
where we see that the reaction density does not depend on $n_{c}$ expressed in
terms of the scale factor. Physically, this can be realized from the fact that
DM production from inflaton scattering depends on the inflaton energy density
and not on the thermal distribution of the bath particles. Finally, we note
that, as this channel is only open during reheating, a boost factor cannot be
defined.

###  5.3 DM from Direct Inflaton Decays

We finally take up the case where the inflaton decay gives rise to the
observed DM abundance through a small branching fraction Br. In this case, the
interaction rate density is

$\gamma=2\,\text{Br}\,\Gamma_{\phi}\,n_{\phi}\simeq
3\,(5-2x)\,\text{Br}\,\frac{M_{P}^{2}\,H_{\text{rh}}^{3}}{m_{\phi}}\left(\frac{T}{T_{\text{rh}}}\right)^{\frac{8\,(3+x)}{3+2x}},$
(5.7)

where $\text{Br}\ll 1$ is required to match the observed DM density, and
therefore the formalism developed in Sec. 3 (i.e., ignoring possible direct
decays of the inflaton into DM) is not modified. With this we obtain the DM
yield at the end of reheating as

$Y(a_{\text{rh}})\simeq\begin{dcases}\frac{3}{2}\,\frac{5-2x}{2x-3}\,\frac{g_{\star}}{g_{\star
s}}\,\text{Br}\,\frac{T_{\text{rh}}}{m_{\phi}}\left[\left(\frac{H_{I}}{H_{\text{rh}}}\right)^{\frac{2x-3}{3}}-1\right]&\text{
for }x\neq\frac{3}{2}\,,\\\ \frac{g_{\star}}{g_{\star
s}}\,\text{Br}\,\frac{T_{\text{rh}}}{m_{\phi}}\ln\left(\frac{H_{I}}{H_{\text{rh}}}\right)&\text{
for }x=\frac{3}{2}\,.\end{dcases}$ (5.8)

Taking into account that
$Y_{0}^{\text{decay}}=\frac{3}{2}\,\frac{g_{\star}}{g_{\star
s}}\,\text{Br}\,\frac{T_{\text{rh}}}{m_{\phi}}$ is the DM yield in the
instantaneous reheating approximation (see, e.g. Ref. [47]), the boost in the
DM yield reads

$B\equiv\frac{Y(a_{\text{rh}})}{Y^{\text{decay}}_{0}}\simeq\begin{dcases}\frac{5-2x}{3-2x}&\text{
for }x<\frac{3}{2}\,,\\\
\frac{2}{3}\ln\left(\frac{H_{I}}{H_{\text{rh}}}\right)&\text{ for
}x=\frac{3}{2}\,,\\\
\frac{5-2x}{2x-3}\left(\frac{H_{I}}{H_{\text{rh}}}\right)^{\frac{2x-3}{3}}&\text{
for }x>\frac{3}{2}\,.\end{dcases}$ (5.9)

Once again, for $x>-3/2$ (or equivalently $k<3/2$), one can write the boost
factor in terms of the temperature ratio $T_{\text{max}}/T_{\text{rh}}$ as

$B\simeq\begin{dcases}\frac{5-2x}{3-2x}&\text{ for }x<\frac{3}{2}\,,\\\
\frac{8}{3+2x}\ln\left(\frac{T_{\text{max}}}{T_{\text{rh}}}\right)&\text{ for
}x=\frac{3}{2}\,,\\\
\frac{5-2x}{2x-3}\left(\frac{T_{\text{max}}}{T_{\text{rh}}}\right)^{\frac{4(2x-3)}{3+2x}}&\text{
for }x>\frac{3}{2}\,.\end{dcases}$ (5.10)

We note that the boost factor can be sizable, i.e., a power law in
$H_{I}/H_{\text{rh}}$ or equivalently $T_{\text{max}}/T_{\text{rh}}$, if
$x>3/2$. On top of that, we also need to have $x<5/2$ for an inflaton
dominated epoch during reheating, as explained before. Hence, to have a
significantly large boost, a large $H_{I}/H_{\text{rh}}$ or equivalently
$T_{\text{max}}/T_{\text{rh}}$ is required for $3/2<x<5/2$, as shown
respectively in the left and right panels of Fig. 4. This figure was produced
using the analytical approximation in Eq. (4.6), in good agreement with the
complete numerical result.

Figure 4: Contours for the boost factor for the DM production from direct
decay of the inflaton, where $B=10^{0}$, $10^{1}$, $10^{2}$, and $10^{3}$ from
bottom left to top right. The gray-shaded area corresponding to $x>5/2$ cannot
be explored in the case of a negligible DM population at the end of inflation.
The left and right panel correspond to the ratios $H_{I}/H_{\text{rh}}$ and
$T_{\text{max}}/T_{\text{rh}}$, respectively.

##  6 Conclusions

In conventional reheating scenarios, a perturbative decay of the inflaton is
typically assumed, which ends at the so-called reheating temperature, whence
the inflaton and the SM radiation have comparable energy densities. In
standard lore, the temperature $T$ of the thermal bath decreases during
reheating as $T(a)\propto a^{-3/8}$ (where $a$ corresponds to the scale
factor). However, this may not have been the case for several instances in
which the temperature of the Universe can either increase or remain steady
during the reheating epoch. For example, it is shown that when the inflaton
oscillates around a potential steeper than quadratic or decays through higher-
order operators, its decay width can feature a time dependence [18, 19, 20,
21, 22, 23], leading to a non-standard evolution of the SM temperature. A
general parametrization of the inflaton dissipation rate can indeed capture
the nontrivial period of reheating. In such a case, the time-dependent decay
of the inflaton not only influences the temperature evolution of the Universe,
but can also critically affect the dark matter (DM) production during
reheating.

By parameterizing the inflaton decay in such a manner, in this work, we have
shown that the UV freeze-in yield of DM can have a power-law boost in the
ratio of the highest temperature and the reheating temperature, namely
$T_{\text{max}}/T_{\text{rh}}$, which is particularly important if the
temperature of the thermal bath drops fast enough during reheating. We find
the evolution of DM yield is controlled by only three free parameters: $x$,
that determines the decay rate of the inflaton, along with the two Hubble
scales $H_{I}$ and $H_{\text{rh}}$ (or equivalently two temperature scales
$T_{\text{max}}$ and $T_{\text{rh}}$). Our main results are derived from a
completely model-agnostic perspective, however we also provide a few instances
where the temperature-dependent inflaton decay can play a nontrivial role in
determining the final DM abundance in context with physically realizable
scenarios. We examine the example of DM production from 2-to-2 scattering of
the inflaton and SM bath particles, where in the former case the yield remains
independent of the inflaton decay dynamics. For DM production through inflaton
decay, we find that a boosted DM yield requires a significantly large
$H_{I}/H_{\text{rh}}$ or equivalently $T_{\text{max}}/T_{\text{rh}}$ ratio
compared to the SM scattering case.

## Acknowledgments

BB and NB received funding from the Patrimonio Autónomo - Fondo Nacional de
Financiamiento para la Ciencia, la Tecnología y la Innovación Francisco José
de Caldas (MinCiencias - Colombia) grant 80740-465-2020. NB and OZ received
funding from MinCiencias through the Grant 80740-492-2021. NB is also funded
by the Spanish FEDER/MCIU-AEI under grant FPA2017-84543-P. The work of OZ has
been partially supported by Sostenibilidad-UdeA, the UdeA/CODI Grants
2017-16286 and 2020-33177. This project has received funding and support from
the European Union’s Horizon 2020 research and innovation programme under the
Marie Skłodowska-Curie grant agreement No. 860881 (H2020-MSCA-ITN-2019
HIDDeN).

## References

  * [1] Planck collaboration, _Planck 2018 results. VI. Cosmological parameters_ , _Astron. Astrophys._ 641 (2020) A6 [1807.06209].
  * [2] G. Bertone and D. Hooper, _History of dark matter_ , _Rev. Mod. Phys._ 90 (2018) 045002 [1605.04909].
  * [3] J. de Swart, G. Bertone and J. van Dongen, _How Dark Matter Came to Matter_ , _Nature Astron._ 1 (2017) 0059 [1703.00013].
  * [4] G. Steigman and M.S. Turner, _Cosmological Constraints on the Properties of Weakly Interacting Massive Particles_ , _Nucl. Phys. B_ 253 (1985) 375.
  * [5] G. Jungman, M. Kamionkowski and K. Griest, _Supersymmetric dark matter_ , _Phys. Rept._ 267 (1996) 195 [hep-ph/9506380].
  * [6] G. Bertone, D. Hooper and J. Silk, _Particle dark matter: Evidence, candidates and constraints_ , _Phys. Rept._ 405 (2005) 279 [hep-ph/0404175].
  * [7] J.L. Feng, _Dark Matter Candidates from Particle Physics and Methods of Detection_ , _Ann. Rev. Astron. Astrophys._ 48 (2010) 495 [1003.0904].
  * [8] L. Roszkowski, E.M. Sessolo and S. Trojanowski, _WIMP dark matter candidates and searches—current status and future prospects_ , _Rept. Prog. Phys._ 81 (2018) 066201 [1707.06277].
  * [9] G. Arcadi, M. Dutra, P. Ghosh, M. Lindner, Y. Mambrini, M. Pierre et al., _The waning of the WIMP? A review of models, searches, and constraints_ , _Eur. Phys. J._ C78 (2018) 203 [1703.07364].
  * [10] J. McDonald, _Thermally generated gauge singlet scalars as selfinteracting dark matter_ , _Phys.Rev.Lett._ 88 (2002) 091304 [hep-ph/0106249].
  * [11] K.-Y. Choi and L. Roszkowski, _E-WIMPs_ , _AIP Conf. Proc._ 805 (2006) 30 [hep-ph/0511003].
  * [12] A. Kusenko, _Sterile neutrinos, dark matter, and the pulsar velocities in models with a Higgs singlet_ , _Phys. Rev. Lett._ 97 (2006) 241301 [hep-ph/0609081].
  * [13] K. Petraki and A. Kusenko, _Dark-matter sterile neutrinos in models with a gauge singlet in the Higgs sector_ , _Phys. Rev._ D77 (2008) 065014 [0711.4646].
  * [14] L.J. Hall, K. Jedamzik, J. March-Russell and S.M. West, _Freeze-In Production of FIMP Dark Matter_ , _JHEP_ 03 (2010) 080 [0911.1120].
  * [15] N. Bernal, M. Heikinheimo, T. Tenkanen, K. Tuominen and V. Vaskonen, _The Dawn of FIMP Dark Matter: A Review of Models and Constraints_ , _Int. J. Mod. Phys._ A32 (2017) 1730023 [1706.07442].
  * [16] F. Elahi, C. Kolda and J. Unwin, _UltraViolet Freeze-in_ , _JHEP_ 03 (2015) 048 [1410.6157].
  * [17] G.F. Giudice, E.W. Kolb and A. Riotto, _Largest temperature of the radiation era and its cosmological implications_ , _Phys. Rev._ D64 (2001) 023508 [hep-ph/0005123].
  * [18] D. Bodeker, _Moduli decay in the hot early Universe_ , _JCAP_ 06 (2006) 027 [hep-ph/0605030].
  * [19] K. Mukaida and K. Nakayama, _Dynamics of oscillating scalar field in thermal environment_ , _JCAP_ 01 (2013) 017 [1208.3399].
  * [20] R.T. Co, E. Gonzalez and K. Harigaya, _Increasing Temperature toward the Completion of Reheating_ , _JCAP_ 11 (2020) 038 [2007.04328].
  * [21] M.A.G. Garcia, K. Kaneta, Y. Mambrini and K.A. Olive, _Inflaton Oscillations and Post-Inflationary Reheating_ , _JCAP_ 04 (2021) 012 [2012.10756].
  * [22] A. Ahmed, B. Grządkowski and A. Socha, _Implications of time-dependent inflaton decay on reheating and dark matter production_ , 2111.06065.
  * [23] A. Banerjee and D. Chowdhury, _Fingerprints of freeze-in dark matter in an early matter-dominated era_ , 2204.03670.
  * [24] T. Moroi and W. Yin, _Light Dark Matter from Inflaton Decay_ , _JHEP_ 03 (2021) 301 [2011.09475].
  * [25] R. Daido, F. Takahashi and W. Yin, _The ALP miracle: unified inflaton and dark matter_ , _JCAP_ 05 (2017) 044 [1702.03284].
  * [26] K. Mukaida, _Dynamics of Scalar Condensation in the Early Universe_ , Ph.D. thesis, Tokyo U., 2015.
  * [27] M.A.G. Garcia, Y. Mambrini, K.A. Olive and M. Peloso, _Enhancement of the Dark Matter Abundance Before Reheating: Applications to Gravitino Dark Matter_ , _Phys. Rev._ D96 (2017) 103510 [1709.01549].
  * [28] N. Bernal, _Boosting Freeze-in through Thermalization_ , _JCAP_ 10 (2020) 006 [2005.08988].
  * [29] N. Bernal, F. Elahi, C. Maldonado and J. Unwin, _Ultraviolet Freeze-in and Non-Standard Cosmologies_ , _JCAP_ 11 (2019) 026 [1909.07992].
  * [30] M.A.G. Garcia, K. Kaneta, Y. Mambrini and K.A. Olive, _Reheating and Post-inflationary Production of Dark Matter_ , _Phys. Rev. D_ 101 (2020) 123507 [2004.08404].
  * [31] N. Bernal, J. Rubio and H. Veermäe, _Boosting Ultraviolet Freeze-in in NO Models_ , _JCAP_ 06 (2020) 047 [2004.13706].
  * [32] R. Allahverdi et al., _The First Three Seconds: a Review of Possible Expansion Histories of the Early Universe_ , _Open J.Astrophys._ 4 (2020)  [2006.16182].
  * [33] S.-L. Chen and Z. Kang, _On UltraViolet Freeze-in Dark Matter during Reheating_ , _JCAP_ 1805 (2018) 036 [1711.02556].
  * [34] N. Bernal, M. Dutra, Y. Mambrini, K. Olive, M. Peloso and M. Pierre, _Spin-2 Portal Dark Matter_ , _Phys. Rev._ D97 (2018) 115020 [1803.01866].
  * [35] G. Bhattacharyya, M. Dutra, Y. Mambrini and M. Pierre, _Freezing-in dark matter through a heavy invisible Z’_ , _Phys. Rev._ D98 (2018) 035038 [1806.00016].
  * [36] D. Chowdhury, E. Dudas, M. Dutra and Y. Mambrini, _Moduli Portal Dark Matter_ , _Phys. Rev._ D99 (2019) 095028 [1811.01947].
  * [37] K. Kaneta, Y. Mambrini and K.A. Olive, _Radiative production of nonthermal dark matter_ , _Phys. Rev._ D99 (2019) 063508 [1901.04449].
  * [38] A. Banerjee, G. Bhattacharyya, D. Chowdhury and Y. Mambrini, _Dark matter seeping through dynamic gauge kinetic mixing_ , _JCAP_ 12 (2019) 009 [1905.11407].
  * [39] P. Chanda, S. Hamdan and J. Unwin, _Reviving $Z$ and Higgs Mediated Dark Matter Models in Matter Dominated Freeze-out_, _JCAP_ 01 (2020) 034 [1911.02616].
  * [40] M. Dutra, _Freeze-in production of dark matter through spin-1 and spin-2 portals_ , _PoS_ LeptonPhoton2019 (2019) 076 [1911.11844].
  * [41] M. Dutra, _The moduli portal to dark matter particles_ , in _11 th International Symposium on Quantum Theory and Symmetries_, 11, 2019, DOI [1911.11862].
  * [42] D. Mahanta and D. Borah, _TeV Scale Leptogenesis with Dark Matter in Non-standard Cosmology_ , _JCAP_ 04 (2020) 032 [1912.09726].
  * [43] C. Cosme, M.a. Dutra, T. Ma, Y. Wu and L. Yang, _Neutrino Portal to FIMP Dark Matter with an Early Matter Era_ , _JHEP_ 03 (2021) 026 [2003.01723].
  * [44] N. Bernal, A. Donini, M.G. Folgado and N. Rius, _Kaluza-Klein FIMP Dark Matter in Warped Extra-Dimensions_ , _JHEP_ 09 (2020) 142 [2004.14403].
  * [45] N. Bernal, J. Rubio and H. Veermäe, _UV Freeze-in in Starobinsky Inflation_ , _JCAP_ 10 (2020) 021 [2006.02442].
  * [46] N. Bernal, A. Donini, M.G. Folgado and N. Rius, _FIMP Dark Matter in Clockwork/Linear Dilaton Extra-Dimensions_ , _JHEP_ 04 (2021) 061 [2012.10453].
  * [47] N. Bernal and Y. Xu, _Polynomial inflation and dark matter_ , _Eur. Phys. J. C_ 81 (2021) 877 [2106.03950].
  * [48] M. Drees, F. Hajkarim and E.R. Schmitz, _The Effects of QCD Equation of State on the Relic Density of WIMP Dark Matter_ , _JCAP_ 1506 (2015) 025 [1503.03513].
  * [49] L. Calibbi, F. D’Eramo, S. Junius, L. Lopez-Honorez and A. Mariotti, _Displaced new physics at colliders and the early universe before its first second_ , _JHEP_ 05 (2021) 234 [2102.06221].
  * [50] R.J. Scherrer and M.S. Turner, _Decaying Particles Do Not Heat Up the Universe_ , _Phys. Rev._ D31 (1985) 681.
  * [51] L. Kofman, A.D. Linde and A.A. Starobinsky, _Towards the theory of reheating after inflation_ , _Phys. Rev. D_ 56 (1997) 3258 [hep-ph/9704452].
  * [52] Planck collaboration, _Planck 2018 results. X. Constraints on inflation_ , _Astron. Astrophys._ 641 (2020) A10 [1807.06211].
  * [53] S. Hannestad, _What is the lowest possible reheating temperature?_ , _Phys. Rev._ D70 (2004) 043506 [astro-ph/0403291].
  * [54] P.F. de Salas, M. Lattanzi, G. Mangano, G. Miele, S. Pastor and O. Pisanti, _Bounds on very low reheating scenarios after Planck_ , _Phys. Rev._ D92 (2015) 123534 [1511.00672].
  * [55] T. Hasegawa, N. Hiroshima, K. Kohri, R.S.L. Hansen, T. Tram and S. Hannestad, _MeV-scale reheating temperature and thermalization of oscillating neutrinos by radiative and hadronic decays of massive particles_ , _JCAP_ 12 (2019) 012 [1908.10189].
  * [56] K. Harigaya, M. Kawasaki, K. Mukaida and M. Yamada, _Dark Matter Production in Late Time Reheating_ , _Phys. Rev._ D89 (2014) 083532 [1402.2846].
  * [57] K. Mukaida and M. Yamada, _Thermalization Process after Inflation and Effective Potential of Scalar Field_ , _JCAP_ 1602 (2016) 003 [1506.07661].
  * [58] K. Harigaya, K. Mukaida and M. Yamada, _Dark Matter Production during the Thermalization Era_ , _JHEP_ 07 (2019) 059 [1901.11027].
  * [59] J. Ellis, M.A.G. Garcia, D.V. Nanopoulos, K.A. Olive and M. Peloso, _Post-Inflationary Gravitino Production Revisited_ , _JCAP_ 1603 (2016) 008 [1512.05701].
  * [60] K. Harigaya and K. Mukaida, _Thermalization after/during Reheating_ , _JHEP_ 05 (2014) 006 [1312.3097].
  * [61] M.A.G. Garcia and M.A. Amin, _Prethermalization production of dark matter_ , _Phys. Rev._ D98 (2018) 103504 [1806.01865].
  * [62] M. Drees and B. Najjari, _Energy spectrum of thermalizing high energy decay products in the early universe_ , _JCAP_ 10 (2021) 009 [2105.01935].
  * [63] M.A.G. Garcia, K. Kaneta, Y. Mambrini, K.A. Olive and S. Verner, _Freeze-in from preheating_ , _JCAP_ 03 (2022) 016 [2109.13280].
  * [64] J.F. Donoghue, _General relativity as an effective field theory: The leading quantum corrections_ , _Phys. Rev. D_ 50 (1994) 3874 [gr-qc/9405057].
  * [65] S.Y. Choi, J.S. Shim and H.S. Song, _Factorization and polarization in linearized gravity_ , _Phys. Rev. D_ 51 (1995) 2751 [hep-th/9411092].
  * [66] B.R. Holstein, _Graviton Physics_ , _Am. J. Phys._ 74 (2006) 1002 [gr-qc/0607045].
  * [67] M. Garny, M. Sandora and M.S. Sloth, _Planckian Interacting Massive Particles as Dark Matter_ , _Phys. Rev. Lett._ 116 (2016) 101302 [1511.03278].
  * [68] Y. Tang and Y.-L. Wu, _Pure Gravitational Dark Matter, Its Mass and Signatures_ , _Phys. Lett. B_ 758 (2016) 402 [1604.04701].
  * [69] M. Garny, A. Palessandro, M. Sandora and M.S. Sloth, _Theory and Phenomenology of Planckian Interacting Massive Particles as Dark Matter_ , _JCAP_ 1802 (2018) 027 [1709.09688].
  * [70] Y. Tang and Y.-L. Wu, _On Thermal Gravitational Contribution to Particle Production and Dark Matter_ , _Phys. Lett. B_ 774 (2017) 676 [1708.05138].
  * [71] Y. Mambrini and K.A. Olive, _Gravitational Production of Dark Matter during Reheating_ , _Phys. Rev. D_ 103 (2021) 115009 [2102.06214].
  * [72] N. Bernal and C.S. Fong, _Dark matter and leptogenesis from gravitational production_ , _JCAP_ 06 (2021) 028 [2103.06896].
  * [73] B. Barman and N. Bernal, _Gravitational SIMPs_ , _JCAP_ 06 (2021) 011 [2104.10699].
  * [74] N. Bernal and Ó. Zapata, _Gravitational dark matter production: primordial black holes and UV freeze-in_ , _Phys. Lett. B_ 815 (2021) 136129 [2011.02510].
  * [75] N. Bernal, _Gravitational Dark Matter and Primordial Black Holes_ , in _Beyond Standard Model: From Theory to Experiment_ , 5, 2021 [2105.04372].
  * [76] A. Ahmed and S. Najjari, _Ultraviolet Freeze-in Dark Matter through the Dilaton Portal_ , 2112.14261.
  * [77] M.R. Haque and D. Maity, _Gravitational dark matter: free streaming and phase space distribution_ , 2112.14668.
  * [78] S. Clery, Y. Mambrini, K.A. Olive and S. Verner, _Gravitational portals in the early Universe_ , _Phys. Rev. D_ 105 (2022) 075005 [2112.15214].
  * [79] M.R. Haque and D. Maity, _Gravitational Reheating_ , 2201.02348.


# Enabling Differentiated Services Using Generalized Power Control Model in
Optical Networks

Quanyan Zhu  Department of Electrical and Computer Engineering
University of Toronto, Ontario M5S 3L1
Email<EMAIL_ADDRESS>Lacra Pavel  Department of Electrical and
Computer Engineering
University of Toronto, Ontario M5S 3L1
Email<EMAIL_ADDRESS>Quanyan Zhu, Lacra Pavel Quanyan Zhu is with
the Department of Electrical and Computer Engineering, University of Illinois
at Urbana Champaign, IL, 61801, USA email<EMAIL_ADDRESS>L. Pavel is
with the Department of Electrical and Computer Engineering, University of
Toronto, Toronto, ON, M5S 3L1 Canada<EMAIL_ADDRESS>
###### Abstract

This paper considers a generalized framework to study OSNR optimization-based
end-to-end link level power control problems in optical networks. We combine
favorable features of game-theoretical approach and central cost approach to
allow different service groups within the network. We develop solutions
concepts for both cases of empty and nonempty feasible sets. In addition, we
derive and prove the convergence of a distributed iterative algorithm for
different classes of users. In the end, we use numerical examples to
illustrate the novel framework.

##  I Introduction

Reconfigurable optical Wavelength-Division Multiplexing (WDM) communication
networks with arbitrary topologies are currently enabled by technological
advances in optical devices such as optical add/drop MUXes (OADM), optical
cross connects (OXC) and dynamic gain equalizer (DGE). It is important that
channel transmission performance and quality of service (QoS) be optimized and
maintained after reconfiguration. At the physical transmission level, channel
performance and QoS are directly determined by the bit-error rate (BER), which
in turn depends on optical signal-to-noise ratio (OSNR), dispersion and
nonlinear effects, [1]. Thus, OSNR is considered as the dominant performance
parameter in link-level optimization. Conventional off-line OSNR optimization
is done by adjusting channel input power at transmitter (Tx) to equalize the
dominant impairment of noise accumulation in chains of optical amplifiers.
However, for reconfigurable optical networks, where different channels can
travel via different optical paths, it is more desirable to implement on-line
decentralized iterative algorithms to accomplish such adjustment.

Recently, this problem is addressed in many research works [2],[3],[4], and
two optimization-based approaches are prevalently used: the central cost and
the non-cooperative game approach. The goals and models of the two approaches
are inherently different. Central cost approach satisfies the target OSNR with
minimum total power consumption. The model embeds the OSNR requirements in its
constraints and indirectly optimizes a certain design criterion. Such model
yields a relatively simple closed-form solution; however, it doesn’t optimize
OSNR in a direct fashion, and thus, channel performance can be potentially
improved for users who need higher quality of transmission. On the other hand,
the game approach is a naturally distributed model which directly optimizes
OSNR based on a payoff function in a non-cooperative manner. Each user
optimizes her own utility to achieve the best possible OSNR. The solution from
this approach is given by Nash equilibrium. As a result, this solution concept
yields best achievable OSNR levels for each user. Since the game approach
involves a cost function arising from pricing, it gives an over-allocation of
resources. Some users may wish to avoid such cost and only demand a basic
level of transmission. Apparently, these two approaches are for two different
type of users and different transmission purposes.

To make use of the advantages from each approach, we propose a generalized
model that combines their features. Such a generalization allows to
accommodate different types of users and also provides a novel mixed framework
to study OSNR power control problem. We separate users into two different
categories. One type of users are those who are willing to pay a price to
fully optimize their transmission performance. Another type of users are those
who are content with basic transmission quality, or OSNR level, set by the
network. The quality of service (QoS) can be met for the former by a game-
theoretically based optimization approach; and for the later by a mechanism
similar to central cost approach.

The contribution of this paper lies in the capability of service
differentiation of the generalized model. For simplicity, total capacity
constraints are not considered. The paper is organized as follows. In section
2, we review the network OSNR model and the basic concepts about the two
optimization-based approaches. In section 3, we establish a general framework
and propose two solution concepts for two different cases of feasible sets.
Section 4 gives an iterative algorithm to achieve such solutions in the
framework. This is illustrated in section 5 by numerical examples. Section 6
concludes the paper and points out future directions of research.

##  II Background

###  II-A Review of Optical Network Model

Consider a network with a set of optical links $\mathcal{L}=\\{1,2,..,L\\}$
connecting the optical nodes, where channel add/drop is realized. A set
$\mathcal{N}=\\{1,2,...,N\\}$ of channels are transmitted, corresponding to a
set of multiplexed wavelengths. Illustrated in Figure 1, a link $l$ has
$K_{l}$ cascaded optically amplified spans. Let $N_{l}$ be the set of channels
transmitted over link $l$. For a channel $i\in\mathcal{N}$, we denote by
$\mathcal{R}_{i}$ its optical path, or collection of links, from source (Tx)
to destination (Rx). Let $u_{i}$ be the $i$th channel input optical power (at
Tx), and $\textbf{u}=[u_{1},...,u_{N}]^{T}$ the vector of all channels’ input
powers. Let $s_{i}$ be the $i$th channel output power (at Rx), and $n_{i}$ the
optical noise power in the $i$th channel bandwidth at Rx. The $i$th channel
optical OSNR is defined as $OSNR_{i}=\frac{s_{i}}{n_{i}}$. In [2], some
assumptions are made to simplify the expression for OSNR, typically for
uniformly designed optical links:

  1. 1.

(A1) Amplified spontaneous emission (ASE) noise power does not participate in
amplifier gain saturation.

  2. 2.

(A2) All the amplifiers in a link have the same spectral shape with the same
total power target and are operated in automatic power control mode.

Under A1 and A2, dispersion and nonlinearity are considered to be limited, and
ASE noise accumulation will be the dominant impairment. The OSNR for the $i$th
channel is given as

$OSNR_{i}=\frac{u_{i}}{n_{0,i}+\sum_{j\in\mathcal{N}}\Gamma_{i,j}u_{j}},i\in\mathcal{N}$
(1)

where $\mathbf{\Gamma}$ is the full $n\times n$ system matrix which
characterizes the coupling between channels. $n_{0,i}$ denotes the $i$th
channel noise power at the transmitter. System matrix $\mathbf{\Gamma}$
encapsulates the basic physics present in optical fiber transmission and
implements an abstraction from a network to an input-output system. This
approach has been used in [5] for the wireless case to model CDMA uplink
communication. Different from the system matrix used in wireless case, the
matrix $\mathbf{\Gamma}$ given in (2) is commonly asymmetric and is more
complicatedly dependent on parameters such as spontaneous emission noise,
wavelength-dependent gain, and the path channels take.

$\Gamma_{i,j}=\sum_{i\in\mathcal{R}_{i}}\sum_{k=1}^{K_{l}}\frac{G_{l,j}^{k}}{G_{l,i}^{k}}\left(\prod_{q=1}^{l-1}\frac{\mathbf{T}_{q,j}}{\mathbf{T}_{q,i}}\right)\frac{ASE_{l,k,i}}{P_{o,l}},\forall
j\in\mathcal{N}_{l}.$ (2)

where $G_{l,k,i}$ is the wavelength dependent gain at $k$th span in $l$th link
for channel $i$; $\mathbf{T}_{l,i}=\prod_{q=1}^{K_{l}}G_{l,k,i}L_{l,k}$ with
$L_{l,k}$ being the wavelength independent loss at $k$th span in $l$th link;
$ASE_{l,k,i}$ is the wavelength dependent spontaneous emission noise;
$P_{0,l}$ is the output power at each span.

Figure 1: A Typical Optical Link in DWDM Optical Networks

###  II-B Central Cost Approach

Similar to the SIR optimization problem in the wireless communication networks
[6, 7], OSNR optimization achieves the target OSNR predefined by each channel
user by allowing the minimum transmission power. Let
$\gamma_{i},i\in\mathcal{N}$ be the target OSNR for each channel. By setting
the OSNR requirement as a constraint, we can arrive at the following central
cost optimization problem (CCP):

(CCP) | $\min_{\mathbf{u}}\sum_{i\in\mathcal{N}}u_{i}$
---|---
subject to | $OSNR_{i}\geq\gamma_{i}\texttt{ }\forall i\in\mathcal{N}.$
(3)

Under certain conditions, it has been shown in [2] that the feasible set of
(CCP) is nonempty and the optimal solution is achievable at the boundary of
the feasible set.

The formulated optimization problem can be extended to incorporate more
constraints such as

$u_{i,\min}\leq u_{i}\leq u_{i,\max},$ (4)

where $u_{i,\min}$ is minimum threshold power required for transmission for
channel $i$ and $u_{i,\max}$ is maximum power channel $i$ can attain. In the
central cost approach, power $u_{i}$ are the parameters to be minimized and
the objective function is linearly separable. In addition, the constraints are
linearly coupled. These nice characteristics in central cost approach leads to
a relatively simple optimization problem.

###  II-C Non-cooperative Game Approach

Let’s review the basic game-theoretical model for power control in optical
networks without constraints. Consider a game defined by a triplet
$\langle\mathcal{N},(A_{i}),(J_{i})\rangle$. $\mathcal{N}$ is the index set of
players or channels; $A_{i}$ is the strategy set $\\{u_{i}\mid
u_{i}\in[u_{i,\min},u_{i,\max}]\\}$; and, $J_{i}$ is the cost function. It is
chosen in a way that minimizing the cost is related to maximizing OSNR level.
In [3], $J_{i}$ is defined as

$J_{i}(u_{i},u_{-i})=\alpha_{i}u_{i}-\beta_{i}\ln\left(1+a_{i}\frac{u_{i}}{X_{-i}}\right),i\in\mathcal{N}$
(5)

where $\alpha_{i},\beta_{i}$ are channel specific parameters, that quantify
the willingness to pay the price and the desire to maximize its OSNR,
respectively, $a_{i}$ is a channel specific parameter, $X_{-i}$ is defined as
$X_{-i}=\sum_{j\neq i}\Gamma_{i,j}u_{j}+n_{0,i}$. This specific choice of
utility function is non-separable, nonlinear and coupled. However, $J_{i}$ is
strictly convex in $u_{i}$ and takes a specially designed form such that its
first-order derivative is linear with respect to $\mathbf{u}$.

The solution from the game approach is usually characterized by Nash
equilibrium (NE). Provided that $\sum_{j\neq i}\Gamma_{i,j}<a_{i}$, the
resulting NE solution is uniquely determined in a closed form by

$\mathbf{\widetilde{{\Gamma}}u^{*}=\widetilde{b}},$ (6)

where $\widetilde{\Gamma}_{i,j}=a_{i},$ for $j=i$;
$\widetilde{\Gamma}_{i,j}=\Gamma_{i,j},$ for $j\neq i$ and
$\widetilde{b}=\frac{a_{i}\beta_{i}}{\alpha_{i}}-n_{0,i}$.

Similar to the wireless case [5], we are able to construct iterative
algorithms to achieve the Nash equilibrium. A simple deterministic first order
parallel update algorithm is:

$u_{i}(n+1)=\frac{\beta_{i}}{\alpha_{i}}-\frac{1}{a_{i}}\left(\frac{1}{OSNR_{i}(n)}-\Gamma_{i,i}\right)u_{i}(n).$
(7)

As proved in [3], the algorithm (7) converges to Nash equilibrium
$\mathbf{u}^{*}$ provided that $\frac{1}{a_{i}}\sum_{j\neq
i}\Gamma_{i,j}<1,\forall i$.

##  III Generalized Model

In this section, we consider a game designed to allow service differentiation
by separating users into two groups: one group seeking a minimum OSNR target
and another group participating in a game setting for OSNR optimization. The
minimum OSNR for target seekers is set by the network to ensure the minimum
quality of service. However, the game players can submit their parameters and
optimize their service accordingly, but they have to pay a price set by the
network for unit power consumption. This concept is illustrated in Figure 2.
Let’s denote set $\mathcal{N}_{1}=\\{1,2,...,N_{1}\\}$ as the set of
competitors, i.e. users who wish to compete for an optimal OSNR. Let set
$\mathcal{N}_{2}=\\{N_{1}+1,\cdots,N_{2}\\}$ be the group of users with target
OSNR given by $\gamma_{i},i\in\mathcal{N}_{2}$. Let
$\mathcal{N}=\mathcal{N}_{1}\cup\mathcal{N}_{2}$, $m=|\mathcal{N}_{1}|=N_{1}$,
$n=|\mathcal{N}_{2}|$, $N=|\mathcal{N}|=m+n$ and
$\mathbf{u}=[u_{1},\cdots,u_{N_{1}},u_{N_{1}+1},\cdots,u_{N_{2}}]^{T}$.

Figure 2: Game players and target seekers in the network

For the game-theoretical players, using the cost function given in (5), we can
form a system of equations given by

$a_{i}u_{i}+X_{-i}=\frac{a_{i}\beta_{i}}{\alpha_{i}},\forall
i\in\mathcal{N}_{1}$

and thus, $\widetilde{\mathbf{\Gamma}}\mathbf{u}=\widetilde{\mathbf{b}},$
where $\widetilde{\mathbf{\Gamma}}\in\mathcal{R}^{m\times N}$ and
$\widetilde{\mathbf{b}}\in\mathcal{R}^{m}$ are defined as in (6). Users with
target OSNR shall have $\mathbf{u}$ satisfy $OSNR_{i}\geq\gamma_{i},\forall
i\in\mathcal{N}_{2},$ or equivalently from (1),

$\frac{u_{i}}{\Gamma_{i,i}u_{i}+\sum_{j\neq
i}\Gamma_{i,j}u_{j}+n_{0,i}}\geq\gamma_{i}$

and thus in a matrix form,
$\widehat{\mathbf{\Gamma}}\mathbf{u}\geq\widehat{\mathbf{b}},$ where
$\widehat{\mathbf{b}}=[\gamma_{1}n_{0,1},\cdots,\gamma_{N}n_{0,N}]^{T}\in\mathcal{R}^{n}$,
$\widehat{\mathbf{\Gamma}}\in\mathcal{R}^{n\times N}$ and is given in (8).

$\widehat{\mathbf{\Gamma}}=\left[\begin{array}[]{ccccc}-\gamma_{N_{1}+1}\Gamma_{N_{1}+1,1}&\cdots&1-\gamma_{N_{1}+1}\Gamma_{N_{1}+1,N_{1}+1}&\cdots&-\gamma_{N_{1}+1}\Gamma_{N_{1}+1,N}\\\
\vdots&\ddots&\ddots&\ddots&\vdots\\\
-\gamma_{N-1}\Gamma_{N-1,1}&-\gamma_{N-1}\Gamma_{N-1,2}&\cdots&1-\gamma_{N-1}\Gamma_{N-1,N-1}&-\gamma_{N-1}\Gamma_{N}\\\
-\gamma_{N}\Gamma_{N,1}&-\gamma_{N}\Gamma_{N,2}&\cdots&\cdots&1-\gamma_{N}\Gamma_{N,N}\\\
\end{array}\right].$ (8)

Let
$F_{1}=\\{\mathbf{u}\in\mathcal{R}^{N}\mid\widetilde{\mathbf{\Gamma}}\mathbf{u}=\widetilde{\mathbf{b}}\\}$
and
$F_{2}=\\{\mathbf{u}\in\mathcal{R}^{N}\mid\widehat{\mathbf{\Gamma}}\mathbf{u}\geq\widehat{\mathbf{b}}\\}$.
In summary, we have a problem formulated as in (DS), where we find solutions
that satisfy $F_{1}$ subject to the constraint described by $F_{2}$.

$\begin{array}[]{cc}\textrm{(DS)}&\widetilde{\mathbf{\Gamma}}\mathbf{u}=\widetilde{\mathbf{b}}\\\
\textrm{s.t.}&\widehat{\mathbf{\Gamma}}\mathbf{u}\geq\widehat{\mathbf{b}}\end{array}$
(9)

In the following discussion, we separate (DS) into two cases: (1) $F=F_{1}\cap
F_{2}\neq\emptyset$, (2)$F=F_{1}\cap F_{2}=\emptyset$, which require different
techniques to find appropriate solutions.

###  III-A Non-empty Feasible Set

A non-empty $F$ may give rise to multiple points that solve (DS). We may
impose some design criteria, or, objective function to reformulate DS for
finding an appropriate solution that solves DS and meet the design criteria at
the same time.

We can use the following result to ensure the nonempty feasible set $F$.

###### Theorem III.1

If
$\overline{\mathbf{\Gamma}}=\left[\begin{array}[]{c}\widetilde{\mathbf{\Gamma}}\\\
\widehat{\mathbf{\Gamma}}\\\ \end{array}\right]$ is nonsingular, the feasible
set $F=F_{1}\cap F_{2}$ is non-empty.

###### Proof:

Let $\mu\in\mathcal{R}^{n}_{+}$ a nonnegative vector. Equivalently, we can
express $F_{2}$ into
$F_{2}=\\{\mathbf{u}\in\mathcal{R}^{n}\mid\widehat{\mathbf{\Gamma}}\mathbf{u}=\widehat{\mathbf{b}}+\mu,\textrm{~{}for~{}some~{}}\mu\in\mathcal{R}^{n}_{+}\\}$.
The set $F$ is thus equivalently
$F=\\{\mathbf{u}\in\mathcal{R}^{N}\mid\overline{\mathbf{\Gamma}}\mathbf{u}=\mathbf{\phi},\textrm{~{}for~{}some~{}}\mu\in\mathcal{R}^{n}_{+}\\}$,
where
$\overline{\mathbf{\Gamma}}=\left[\begin{array}[]{c}\widetilde{\mathbf{\Gamma}}\\\
\widehat{\mathbf{\Gamma}}\\\ \end{array}\right]$ and
$\mathbf{\phi}=\left[\begin{array}[]{c}\widetilde{\mathbf{b}}\\\
\widehat{\mathbf{b}}+\mu\\\ \end{array}\right]$. If
$\overline{\mathbf{\Gamma}}$ is nonsingular, there exist a unique
$\mathbf{u}\in\mathcal{R}^{N}$ for every nonnegative $\mu$. Therefore $F$ is
non-empty. ∎

Suppose conditions in Theorem III.1 hold and $F$ is nonempty. We consider an
appropriate solution in $F$ that satisfies a certain design criteria. Thus, we
formulate (DSNP111DSNP stands for “Differentiated Service N-person Problem”.)
in which we minimize total power consumption subject to the conditions arising
from the different service requirements.

$\begin{array}[]{cc}\textrm{(DSNP)}&\min\sum_{i}u_{i}\\\
\textrm{s.t.}&\widetilde{\mathbf{\Gamma}}\mathbf{u}=\widetilde{\mathbf{b}},\widehat{\mathbf{\Gamma}}\mathbf{u}\geq\widehat{\mathbf{b}}\end{array}$
(10)

The constraints of (DSNP) can be relaxed and augmented into

$\overline{\mathbf{\Gamma}}\mathbf{u}\geq\overline{\mathbf{b}}.$ (11)

where
$\overline{\mathbf{\Gamma}}=\left[\begin{array}[]{c}\widetilde{\mathbf{\Gamma}}\\\
\widehat{\mathbf{\Gamma}}\\\ \end{array}\right]\in\mathcal{R}^{N\times N}$ and
$\overline{\mathbf{b}}=\left[\begin{array}[]{c}\widetilde{\mathbf{b}}\\\
\widehat{\mathbf{b}}\\\ \end{array}\right]\in\mathcal{R}^{N}.$

According to the fundamental theorem of linear programming [8], if (DSNP) is
realistic, the solution is obtained at the extreme point of the feasible set
$F$. Since $F$ has only one extreme point when $\overline{\mathbf{\Gamma}}$ is
non-singular, the solution is uniquely given by

$\mathbf{u}=\overline{\mathbf{\Gamma}}^{-1}\overline{\mathbf{b}}.$ (12)

To further characterize the solution $\mathbf{u}$, we assume strict diagonal
dominance of matrix $\overline{\mathbf{\Gamma}}$ [9], which leads to non-
singularity of the matrix and uniqueness of the solution.

###### Theorem III.2

Suppose OSNR targets $\gamma_{i},i\in\mathcal{N}_{2}$ are chosen such that
$\gamma_{i}<\frac{1}{\sum_{j\in\mathcal{N}}\Gamma_{i,j}},i\in\mathcal{N}_{2}$.
In addition, parameters $a_{i}$ are chosen as $a_{i}>\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{ij},\forall i\in\mathcal{N}_{1}.$ The matrix
$\overline{\mathbf{\Gamma}}$ is strictly diagonally dominant. And thus, a
unique solution to (DSNP) is given by (12).

###### Proof:

From the assumption that
$\gamma_{i}\sum_{j\in\mathcal{N}}\Gamma_{ij}<1,i\in\mathcal{N}_{2}$, it is
apparent that $\gamma_{i}<\frac{1}{\Gamma_{ii}}$ and
$\left|1-\gamma_{i}\Gamma_{ii}\right|>\gamma_{i}\sum_{j}\Gamma_{ij},\forall
i\in\mathcal{N}_{2}$. In addition, $a_{i}>\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j},\forall i\in\mathcal{N}_{1}$. Therefore, matrix
$\overline{\mathbf{\Gamma}}$ is strictly diagonally dominant. Using Gershgorin
theorem in [9], we conclude that there exists a unique solution to (DSNP). ∎

The assumption of strict diagonal dominance in Theorem III.2 is reasonable
because typical values of $\Gamma_{ij}$ are found to be on the order of
$10^{-3}$ and desirable levels of OSNR are 20-30dB.

###### Remark III.1

(DSNP) can be seen as a generalized approach that combines central cost
approach in [2] and non-cooperative game approach in [3]. When
$N_{1}=\emptyset,N_{2}\neq\emptyset$, (DSNP) reduces to the central cost
approach. Similarly, when $N_{1}\neq\emptyset,N_{2}=\emptyset$, (DSNP) reduces
to the game-theoretical approach and the given solution is Nash equilibrium
accordingly. This framework allows to study two different types of users at
the same time.

###### Remark III.2

We illustrate a two-person (DSNP), where player 1 chooses to compete and
optimize his utility and player 2 chooses to meet a certain OSNR target
$\gamma_{2}$. We form the 2-by-2 matrix $\overline{\mathbf{\Gamma}}$ and
$\overline{\mathbf{b}}$ as follows.

$\overline{\mathbf{\Gamma}}=\left[\begin{array}[]{cc}a_{1}&\Gamma_{12}\\\
-\Gamma_{21}\gamma_{2}&1-\Gamma_{22}\gamma_{2}\\\
\end{array}\right],\overline{\mathbf{b}}=\left[\begin{array}[]{c}\frac{a_{1}\beta_{1}}{\alpha_{1}}-n_{0,1}\\\
n_{0,2}\gamma_{2}\\\ \end{array}\right]$

The feasible set $F=F_{1}\cap F_{2}$ is shown in Figure 3 by a dotted line.
The relaxed (DSNP) has its relaxed feasible depicted in the shaded region. The
solution is given by
$\mathbf{u}^{*}=\overline{\mathbf{\Gamma}}^{-1}\overline{\mathbf{b}},$ which
is illustrated by the dark point in Figure 3. $\mathbf{u}^{*}$ is nonnegative
componentwise if network price $\alpha_{1}$ is set such that
$s_{2}>\frac{n_{0,2}}{1-\Gamma_{22}}$.

Figure 3: The feasible set of two-person (DSNP).
$s_{1}=\frac{\tilde{b}_{1}}{a_{1}}$; $s_{2}=\frac{\tilde{b}_{1}}{\Gamma_{12}}$

Based on Theorem III.2, we can further investigate how parameters chosen by
game players and target seekers influence the outcome of the allocation. The
result is summarized in Theorem III.3.

###### Theorem III.3

Let $\kappa$ be the condition number of $\overline{\mathbf{\Gamma}}$,
$T_{i}=a_{i}+\sum_{j\neq i,j\in\mathcal{N}}\Gamma_{ij},\forall
i\in\mathcal{N}_{1}$ and $S_{k}=2-2\gamma_{k}\Gamma_{kk},\forall
k\in\mathcal{N}_{2}$. Suppose $\overline{\mathbf{\Gamma}}$ is strictly
diagonally dominant by satisfying conditions in Theorem III.2. In addition,
$T_{i}>S_{k}$ and $\tilde{b}_{i}>\hat{b}_{k},\forall
i\in\mathcal{N}_{1},\forall k\in\mathcal{N}_{2}.$ The maximum allocated power
allocated to users are bound as follows.

$\frac{\max_{i\in\mathcal{N}_{2}}\gamma_{i}n_{0,i}}{\max_{i\in\mathcal{N}_{1}}2a_{i}}\leq\|\mathbf{u}\|_{\infty}\leq\kappa\max_{i\in\mathcal{N}_{1}}\frac{\beta_{i}}{\alpha_{i}}$

###### Proof:

Let $R_{i}$ denote the i-th row absolute sum of matrix
$\overline{\mathbf{\Gamma}}$, i.e.,

$R_{i}=\sum_{j\in\mathcal{N}}\left|\overline{\Gamma}_{ij}\right|.$ (13)

Using conditions from Theorem III.2, we arrive at

$R_{i}=\left\\{\begin{array}[]{ll}1+\gamma_{i}\sum_{j\in\mathcal{N}}\Gamma_{ij}-2\gamma_{i}\Gamma_{ii}<2-2\gamma_{i}\Gamma_{ii},&{i\in\mathcal{N}_{2};}\\\
a_{i}+\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{ij}<2a_{i}.,&{i\in\mathcal{N}_{1}.}\end{array}\right.$
(14)

With the assumption that $a_{i}+\sum_{j\neq
i}\Gamma_{ij}>2-2\gamma_{k}\Gamma_{kk},\forall i\in\mathcal{N}_{1},\forall
k\in\mathcal{N}_{2},$ we obtain
$\|\overline{\mathbf{\Gamma}}\|_{\infty}=\max_{i\in\mathcal{N}}R_{i}=\max_{i\in\mathcal{N}}a_{i}+\sum_{j\neq
i}\Gamma_{ij}.$ Using (14) and the fact that $\Gamma_{ij}\geq 0$, we obtain an
upper and lower bound on $\|\overline{\mathbf{\Gamma}}\|_{\infty}$, i.e.,

$\max_{i\in\mathcal{N}_{1}}a_{i}\leq\|\overline{\mathbf{\Gamma}}\|_{\infty}\leq\max_{i\in\mathcal{N}_{1}}2a_{i}.$
(15)

In addition, from $\tilde{b}_{i}>\hat{b}_{k},\forall
i\in\mathcal{N}_{1},\forall k\in\mathcal{N}_{2},$ we obtain an upper bound and
lower bound for $\|\overline{\mathbf{b}}\|_{\infty}$, given by

$\max_{i\in\mathcal{N}_{2}}\gamma_{i}n_{0,i}\leq\|\overline{\mathbf{b}}\|_{\infty}=\max_{i\in\mathcal{N}}\overline{b}_{i}\leq\max_{i\in\mathcal{N}_{1}}\tilde{b}_{i}=\max_{i\in\mathcal{N}_{1}}\frac{a_{i}\beta_{i}}{\alpha_{i}}$
(16)

Since $\overline{\mathbf{\Gamma}}$ is strictly diagonally dominant, using
matrix norm sub-multiplicativity, we obtain from (12)

$\frac{\|\overline{\mathbf{b}}\|_{\infty}}{\|\overline{\mathbf{\Gamma}}\|_{\infty}}\leq\|\mathbf{u}\|_{\infty}\leq\frac{\kappa\|\overline{\mathbf{b}}\|_{\infty}}{\|\overline{\mathbf{\Gamma}}\|_{\infty}},$
(17)

where $\kappa$ is the condition number of $\overline{\mathbf{\Gamma}}$ given
by
$\kappa=\|\overline{\mathbf{\Gamma}}\|_{\infty}\|\overline{\mathbf{\Gamma}}^{-1}\|_{\infty}\geq
1.$

Using (15), (16) and (17), we obtain

$\displaystyle\frac{\max_{i\in\mathcal{N}_{2}}\gamma_{i}n_{0,i}}{\max_{i\in\mathcal{N}_{1}}2a_{i}}\leq\|\mathbf{u}\|_{\infty}$
$\displaystyle\leq$
$\displaystyle\frac{\kappa\max_{i\in\mathcal{N}_{1}}a_{i}\beta_{i}/\alpha_{i}}{\max_{i\in\mathcal{N}_{1}}a_{i}}$
(18) $\displaystyle\leq$
$\displaystyle\frac{\kappa\max_{i\in\mathcal{N}_{1}}a_{i}\max_{i\in\mathcal{N}_{1}}\beta_{i}/\alpha_{i}}{\max_{i\in\mathcal{N}_{1}}a_{i}}$
$\displaystyle\leq$
$\displaystyle\kappa\max_{i\in\mathcal{N}_{1}}\frac{\beta_{i}}{\alpha_{i}}.$

∎

It is easy to observe that the upper bound is dependent on the parameters of
the game players and the lower bound is dependent on the OSNR levels of target
seeker and parameter $a_{i}$ of the game players. In essence, game players
control the outcome of the model and the choice of OSNR target can only affect
the lower bound. Such relation describes a fair scenario in which game
players, who pay for their power at $\alpha_{i}$, have their choices of
parameters $a_{i},\beta_{i}$ to influence the network allocation.

###### Remark III.3

Since
$\|\mathbf{u}\|_{\infty}\leq\|\mathbf{u}\|_{2}\leq\sqrt{N}\|\mathbf{u}\|_{\infty}$,
we can translate the result obtained in (18) directly into Euclidean norm,
i.e.,

$B_{\infty}^{L}\leq\|\mathbf{u}\|_{2}\leq\sqrt{N}B^{U}_{\infty}$ (19)

where
$B_{\infty}^{U}=\kappa\max_{i\in\mathcal{N}_{1}}\frac{\beta_{i}}{\alpha_{i}}$
and
$B_{\infty}^{L}=\frac{\max_{i\in\mathcal{N}_{2}}\gamma_{i}n_{0,i}}{\max_{i\in\mathcal{N}_{1}}2a_{i}}$.

By (19), we can see that the network can encourage uniform channel power
distribution by letting $B_{\infty}^{U}$ close to $\sqrt{N}B_{\infty}^{L}$ and
provide incentive for differentiated services by letting them far apart. It
can be implemented by the network by adjusting OSNR level $\gamma_{i}$ and
pricing $\alpha_{i}$. Decreasing $\alpha_{i}$ encourages more users to be game
players, giving rise to more competitions or service differentiation as a
result of higher upper bound. On the other hand, increasing $\gamma_{i}$
raises the lower bound and encourages more users being target-seekers.

###  III-B Empty Feasible Set

In this section, we consider the second case where feasible set $F$ is empty.
Instead of finding an appropriate feasible solution, we find the closest
points between set $F_{1}$ and $F_{2}$. We use a quadratic program (DS2) to
minimize the error norm subject to the constraint described by $F_{2}$.

$\begin{array}[]{cc}\textrm{(DS2)}&\min_{\mathbf{u}}\|\widetilde{\mathbf{\Gamma}}\mathbf{u}-\widetilde{\mathbf{b}}\|_{2}\\\
\textrm{s.t.}&\widehat{\mathbf{\Gamma}}\mathbf{u}\geq\widehat{\mathbf{b}}\end{array}$
(20)

We can turn the constrained problem (20) into an unconstrained problem by
studying its corresponding dual problem. Since
$\|\widetilde{\mathbf{\Gamma}}\mathbf{u}-\widetilde{\mathbf{b}}\|_{2}=\mathbf{u}^{T}\widetilde{\mathbf{\Gamma}}^{T}\widetilde{\mathbf{\Gamma}}\mathbf{u}-2(\widetilde{\mathbf{b}}^{T}\widetilde{\mathbf{\Gamma}})\mathbf{u}+\widetilde{\mathbf{b}}^{T}\widetilde{\mathbf{b}}$,
we denote
$\mathbf{H}=\frac{1}{2}\widetilde{\mathbf{\Gamma}}^{T}\widetilde{\mathbf{\Gamma}},\mathbf{d}=-2(\widetilde{\mathbf{\Gamma}}^{T}\widetilde{\mathbf{b}})$,
$\mathbf{D}=-\widehat{\mathbf{\Gamma}}(\mathbf{H}^{T}\mathbf{H})^{-1}\mathbf{H}^{T}\widehat{\mathbf{\Gamma}}^{T}$,
$\mathbf{c}=\widehat{\mathbf{b}}+\widehat{\mathbf{\Gamma}}(\mathbf{H}^{T}\mathbf{H})^{-1}\mathbf{H}^{T}\mathbf{d}$;
and form a Lagrangian from the original problem (DS2).

$\displaystyle D(\mu)$ $\displaystyle=$
$\displaystyle\min_{\mathbf{u}}\mathcal{L}(\mathbf{u},\mu)$ $\displaystyle=$
$\displaystyle\min_{\mathbf{u}}\left(\frac{1}{2}\mathbf{u}^{T}\mathbf{H}\mathbf{u}+\mathbf{d}^{T}\mathbf{u}+\widetilde{\mathbf{b}}^{T}\widetilde{\mathbf{b}}+\mu^{T}(-\widehat{\mathbf{\Gamma}}\mathbf{u}+\widetilde{\mathbf{b}})\right)$

Since the objective function is convex, the necessary and sufficient condition
for a minimum is that the gradient must vanish,i.e.,

$\mathbf{H}\mathbf{u}+\mathbf{d}-\hat{\mathbf{\Gamma}}^{T}\mathbf{\mu}=0.$
(22)

For $n<N$, $\widetilde{\mathbf{\Gamma}}$ is not full rank. Therefore,
$\mathbf{H}$ is singular and there exist multiple solutions to (22). Using
pseudoinverse [9], we can find a solution to (22) given by

$\mathbf{u}=-(\mathbf{H}^{T}\mathbf{H})^{-1}\mathbf{H}^{T}\left(\mathbf{d}-\hat{\mathbf{\Gamma}}^{T}\mu\right).$

Thus, after replacing into (III-B), we obtain $\mu$ as a solution to the dual
problem (DDS2).

$\textrm{(DDS2)}\max_{\mu\geq
0}\frac{1}{2}\mu^{T}\mathbf{D}\mu+\mu^{T}\mathbf{c}-\frac{1}{2}\mathbf{d}^{T}(\mathbf{H}^{T}\mathbf{H})^{-1}\mathbf{H}^{T}\mathbf{d}+\mathbf{b}^{T}\mathbf{b}$
(23)

The problem (LDS2) and dual problem (DDS2) can be solved using unconstrained
optimization algorithms in [10], [8].

##  IV Iterative Algorithm

In this section, we develop algorithm for the case of nonempty $F$ set. Let
$u_{i}(n)$ denote the power at channel $i$ at step $n$. An iterative algorithm
is given as follows.

$\left\\{\begin{array}[]{ll}u_{i}(n+1)=\frac{\beta_{i}}{\alpha_{i}}-\frac{1}{a_{i}}\left(\frac{1}{OSNR_{i}(n)}-\Gamma_{i,i}\right)u_{i}(n),&\forall
i\in\mathcal{N}_{1};\\\
u_{i}(n+1)=\frac{\gamma_{i}}{1-\gamma_{i}\Gamma_{i,i}}\left(\frac{1}{OSNR_{i}(n)}-\Gamma_{i,i}\right)u_{i}(n),&\forall
i\in\mathcal{N}_{2}.\end{array}\right.$ (24)

###### Theorem IV.1

Algorithm (24) converges provided that $a_{i}>\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}$ and $\gamma_{i}$ is chosen such that
$\gamma_{i}<\frac{1}{\sum_{j\in\mathcal{N}}\Gamma_{i,j}}$.

###### Proof:

We use a similar approach from [3] to show the convergence of (24). Let’s
define $e_{i}(n)=u_{i}(n)-u_{i}^{*}$, where $u_{i}^{*}$ is given in (12).
Since $\overline{\mathbf{\Gamma}}\mathbf{u}^{*}=\overline{\mathbf{b}}$,
$\widetilde{\Gamma}_{i,i}u_{i}^{*}+\sum_{j\neq
i}\widetilde{\Gamma}_{i,j}u_{j}^{*}=\tilde{b}_{i}$, for $i\in\mathcal{N}_{1}$;
and, $\widehat{\Gamma}_{i,i}u_{i}^{*}+\sum_{j\neq
i}\widehat{\Gamma}_{i,j}u_{j}^{*}=\hat{b}_{i}$, for $i\in\mathcal{N}_{2}$.

Substitute the expression for $u_{i}^{*}$ into $e_{i}(n+1)$, and we obtain
$e_{i}(n+1)=u_{i}(n+1)-u_{i}^{*}=-\frac{1}{{a_{i}}}\left[\sum_{j\neq
i}{\Gamma}_{i,j}(u_{j}(n)-u_{j}^{*})\right]$, for $i\in\mathcal{N}_{1}$; and
$e_{i}(n+1)=u_{i}(n+1)-u_{i}^{*}=\frac{1}{1-{\Gamma}_{i,i}\gamma_{i}}\left[\sum_{j\neq
i}{\Gamma}_{i,j}\gamma_{i}(u_{j}(n)-u_{j}^{*})\right]$, for
$i\in\mathcal{N}_{2}$. Let $\mathbf{e}=[e_{i}(n)],i\in\mathcal{N}$. Therefore,
for $i\in\mathcal{N}_{1}$,

$\displaystyle|e_{i}(n+1)|$ $\displaystyle=$
$\displaystyle\left|\frac{1}{a_{i}}\left[\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}(e_{j}(n))\right]\right|$ (25)
$\displaystyle\leq$ $\displaystyle\frac{1}{a_{i}}\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}\max_{j\in\mathcal{N}}|e_{j}(n)|$ (26)
$\displaystyle\leq$ $\displaystyle\frac{1}{a_{i}}\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}\|\mathbf{e}(n)\|_{\infty}.$ (27)

and similarly, for $i\in\mathcal{N}_{2}$,

$\displaystyle|e_{i}(n+1)|$ $\displaystyle=$
$\displaystyle\left|\frac{1}{1-\Gamma_{i,i}\gamma_{i}}\left[\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}\gamma_{i}(e_{j}(n))\right]\right|$ (28)
$\displaystyle\leq$
$\displaystyle\frac{\gamma_{i}}{|1-\Gamma_{i,i}\gamma_{i}|}\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}\max_{j\in\mathcal{N}}|e_{j}(n)|.$
$\displaystyle\leq$
$\displaystyle\frac{\gamma_{i}}{|1-\Gamma_{i,i}\gamma_{i}|}\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}\|\mathbf{e}(n)\|_{\infty}.$ (29)

Since we assumed that $a_{i}>\sum_{j\neq i,j\in\mathcal{N}}\Gamma_{i,j}$ and
$\gamma_{i}$ is chosen such that
$\gamma_{i}<\frac{1}{\sum_{j\in\mathcal{N}}\Gamma_{i,j}}\leq\frac{1}{\Gamma_{i,i}}$,
we can conclude that $\|\mathbf{e}(n)\|\rightarrow 0$ from the contraction
mapping theorem. As a result, we have $u_{i}(n)\rightarrow u_{i}^{*}$ as
$n\rightarrow\infty$, for $i\in\mathcal{N}$. ∎

###### Remark IV.1

From the proof, we note that the rate of convergence of is determined by

$\sigma=\max\left\\{\max_{i\in\mathcal{N}_{1}}\frac{\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}}{a_{i}},\max_{i\in\mathcal{N}_{2}}\frac{\sum_{j\neq
i,j\in\mathcal{N}}\Gamma_{i,j}\gamma_{i}}{1-\Gamma_{i,i}\gamma_{i}}\right\\}.$

In addition, it is easy to observe that the OSNR target-seeking users are
algorithmically equivalent to competition seeking users by letting
$\beta_{i}/\alpha_{i}=0$ and $a_{i}=\Gamma_{i,i}-\frac{1}{\gamma_{i}}$,
$i\in\mathcal{N}_{2}$. This is because no notion of pricing is used for the
OSNR target seekers and they just have a utility target to meet or
equivalently optimize by letting $a_{i}=\Gamma_{i,i}-\frac{1}{\gamma_{i}}$.

##  V Numerical Examples

In this section, we illustrate the concept by a MATLAB simulation. We consider
an end-to-end link described in Figure 1 with 5 amplified spans. We assume
channels are transmitted at wavelengths distributed centered around 1555nm
with channel separation of 1nm. Suppose input noise power is 0.5 percent of
the input signal power. The gain profile for each amplifier is identically
assumed to be parabolic as in Figure 4, which is normalized with respect to
$G_{\max}=30.0$dB. Suppose 20dB is the target OSNR level for users who just
want to meet a satisfactory level of transmission. We first show the case of 3
users, in which 2 users need better quality of service and one user is simply
interested in meeting 20dB as a target. From Figure 5, we can observe that
users who need better services reach an OSNR of 26.33dB and 29.20dB,
respectively. With an appropriate choice of initial conditions, the algorithm
quickly converges in 1-2 steps. In Figure 6, we similarly show the case of 30
users, in which 20 are game players and 10 are target seekers.

Figure 4: Optical Amplifier Spectral Profile Figure 5: OSNR simulation with 3
users in time steps Figure 6: OSNR simulation with 30 users in time steps

##  VI Conclusion

In this paper, we examined a generalized power control model in optical
networks, which combines features of central cost approach and game-
theoretical approach. It enables two major service types in the network. One
is game player, who pays for his power consumption and the other is target
seeker, who is satisfied with a minimum service level set by the network. We
discussed two different solutions concepts for nonempty and empty feasible set
respectively and specifically designed an iterative algorithm that converges
to a unique solution for the case of nonempty feasible set. The convergence of
the algorithm was proved and illustrated by numerical examples of a WDM end-
to-end optical link.

In this work, we didn’t include capacity constraints for the sake of
simplicity. We hope this work will lead to future investigations of more
complicated cases where network constraints and nonlinear effects are
considered. In addition, we expect this framework to be used to solve similar
problems in other types of networks, for example, wireless networks.

## References

  * [1] G. Agrawal, _Lightwave Technology_. Wiley-Interscience, 2005.
  * [2] L. Pavel, “OSNR optimization in optical networks: Modeling and distributed algorithms via a central cost approach,” _IEEE Journal on Selected Areas in Communications_ , vol. 24, no. 4, pp. 54–65, April 2006.
  * [3] ——, “A noncooperative game approach to OSNR optimization in optical networks,” _IEEE Transactions on Automatic Control_ , vol. 51, no. 5, pp. 848–852, May 2006.
  * [4] Y. Pan and L. Pavel, “OSNR optimization in optical networks: Extension for capacity constraints,” _Proceedings of 2005 American Control Conference_ , pp. 2379–2385, June 2005.
  * [5] R. Srikant, E. Altman, T. Alpcan, and T. Basar, “CDMA uplink power control as noncooperative game,” _Wireless Networks_ , vol. 8, p. 659 690, 2002\.
  * [6] C. Saraydar, N. Mandayam, and D. Goodman, “Efficient power control via pricing in wireless data networks,” _IEEE Transactions on Communications_ , vol. 50, no. 2, pp. 291–414, February 2002.
  * [7] C. Saraydar and D. Goodman, “Pricing and power control in a multicell wireless data network,” _IEEE Journal of Selected Areas of Communications_ , vol. 19, no. 10.
  * [8] D. Bertsekas, _Nonlinear Programming_. Athena Scientific, 2003.
  * [9] R. Horn and C. Johnson, _Matrix Analysis_. Cambridge University Press, 1990.
  * [10] M. Bazaraa, H. Sherali, and C. Shetty, _Nonlinear Programming: Theory and Algorithms_ , 2nd ed. Wiley, 1993.


* Guastaroba et al. [2011] Guastaroba, G., Mitra, G., and Speranza, M. G. (2011).  Investigating the effectiveness of robust portfolio optimization techniques.  Journal of Asset Management, 12(4):260–280.
  * Kim and Boyd [2008] Kim, S.-J. and Boyd, S. (2008).  A minimax theorem with applications to machine learning, signal processing, and finance.  SIAM Journal on Optimization, 19(3):1344–1367.
  * Lin [1991] Lin, J. (1991).  Divergence measures based on the shannon entropy.  IEEE Transactions on Information theory, 37(1):145–151.
  * Lobo and Boyd [2000] Lobo, M. S. and Boyd, S. (2000).  The worst-case risk of a portfolio.  Technical report. Available from http://web.stanford.edu/$\sim$boyd/papers/pdf/risk_bnd.pdf.
  * Maillard et al. [2010] Maillard, S., Roncalli, T., and Teiletche, J. (2010).  The properties of equally weighted risk contribution portfolios.  Journal of Portfolio Management, 36(4):60–70.
  * Markowitz [1952] Markowitz, H. (1952).  Portfolio selection.  Journal of Finance, 7(1):77–91.
  * Mausser and Romanko [2014] Mausser, H. and Romanko, O. (2014).  Computing equal risk contribution portfolios.  IBM Journal of Research and Development, 58(4):5–1.
  * Merton [1980] Merton, R. C. (1980).  On estimating the expected return on the market: An exploratory investigation.  Journal of financial economics, 8(4):323–361.
  * Nedić and Ozdaglar [2009] Nedić, A. and Ozdaglar, A. (2009).  Subgradient methods for saddle-point problems.  Journal of optimization theory and applications, 142(1):205–228.
  * Rahimian and Mehrotra [2019] Rahimian, H. and Mehrotra, S. (2019).  Distributionally robust optimization: A review.  arXiv preprint arXiv:1908.05659.
  * Rustem and Howe [2009] Rustem, B. and Howe, M. (2009).  Algorithms for worst-case design and applications to risk management.  Princeton University Press.
  * Scarf [1958] Scarf, H. (1958).  A min-max solution of an inventory problem.  Studies in the mathematical theory of inventory and production, pages 201–209.
  * Shapiro and Ahmed [2004] Shapiro, A. and Ahmed, S. (2004).  On a class of minimax stochastic programs.  SIAM Journal on Optimization, 14(4):1237–1249.
  * Shapiro et al. [2014] Shapiro, A., Dentcheva, D., and Ruszczyński, A. (2014).  Lectures on stochastic programming: modeling and theory.  SIAM.
  * Shapiro and Kleywegt [2002] Shapiro, A. and Kleywegt, A. (2002).  Minimax analysis of stochastic problems.  Optimization Methods and Software, 17(3):523–542.
  * Sharpe [1964] Sharpe, W. F. (1964).  Capital asset prices: A theory of market equilibrium under conditions of risk.  The journal of finance, 19(3):425–442.
  * Sharpe [1994] Sharpe, W. F. (1994).  The sharpe ratio.  Journal of Portfolio Management, 21(1):49–58.
  * Tütüncü and Koenig [2004] Tütüncü, R. H. and Koenig, M. (2004).  Robust asset allocation.  Annals of Operations Research, 132(1-4):157–187.
  * Von Neumann [1928] Von Neumann, J. (1928).  Zur theorie der gesellschaftsspiele.  Mathematische annalen, 100(1):295–320.
  * Xiu and Zhang [2003] Xiu, N. and Zhang, J. (2003).  Some recent advances in projection-type methods for variational inequalities.  Journal of Computational and Applied Mathematics, 152(1-2):559–585.
  * Žáčková [1966] Žáčková, J. (1966).  On minimax solutions of stochastic linear programming problems.  Časopis pro pěstování matematiky, 91(4):423–430.

##  Appendix A Numerical implementation of statistical distances

Here we describe how to numerically implement the squared Hellinger distance
in (16) and the TV distance in (17). We use either of these two distance
measures to define the ambiguity set $\mathcal{U}_{\bm{p}}$, and then use the
set to construct the corresponding Euclidean projection optimization problem
in (29). However, in their current form, most optimization solvers will reject
them.

If we wish to use the squared Hellinger distance in (16), then the projection
optimization problem can be implemented as follows.

$\displaystyle\begin{aligned}
&\min_{\bm{p},\bm{r}}&\|\bm{u}-\bm{p}\|_{2}^{2}\end{aligned}$
$\displaystyle\begin{aligned} &\ \text{s.t.}&\bm{1}^{T}\bm{p}&=1,\\\
&&\frac{1}{2}\sum_{t=1}^{T}p_{t}-2r_{t}\sqrt{q_{t}}+q_{t}&\leq
d_{\text{H}},\\\ &&p_{t}&\geq r_{t}^{2}\quad\text{for }t=1,\dots,T,\\\
&&\bm{p},\ \bm{r}&\geq 0,\end{aligned}$

where $\bm{u}\in\mathbb{R}^{T}$ is some arbitrary vector that we wish to
project onto the set $\mathcal{U}_{\bm{p}}$, while $\bm{r}\in\mathbb{R}^{T}$
is an auxiliary variable that serves as a placeholder for the square root of
each element of $\bm{p}$. As before, $\bm{q}\in\mathcal{P}$ is the nominal
probability distribution, while $d_{\text{H}}$ is the maximum permissible
distance in (20) and is defined by the number of scenarios $T$ and the
investor’s desired degree of robustness $\omega$.

On the other hand, if we wish to use the TV distance in (17), then the
projection optimization problem can be implemented as follows.

$\displaystyle\begin{aligned}
&\min_{\bm{p},\bm{\zeta}}&\|\bm{u}-\bm{p}\|_{2}^{2}\end{aligned}$
$\displaystyle\begin{aligned} &\ \text{s.t.}&\bm{1}^{T}\bm{p}&=1,\\\
&&\frac{1}{2}\sum_{t=1}^{T}\zeta_{t}&\leq d_{\text{TV}},\\\ &&\zeta_{t}&\geq
p_{t}-q_{t}\quad\text{for }t=1,\dots,T,\\\ &&\zeta_{t}&\geq
q_{t}-p_{t}\quad\text{for }t=1,\dots,T,\\\ &&\bm{p}&\geq 0,\end{aligned}$

where $\bm{\zeta}\in\mathbb{R}^{T}$ is an auxiliary variable that represents
the absolute value of the difference between the elements of $\bm{p}$ and
$\bm{q}$. The maximum permissible distance $d_{\text{TV}}$ is defined by the
number of scenarios $T$ and the investor’s desired degree of robustness
$\omega$ as shown in (21).