# Structure of hot strange quark stars: an NJL model approach at finite
temperature

G. H. Bordbar1,2,3 111Corresponding author. E-mail<EMAIL_ADDRESS>R. Hosseini 1 and F. Kayanikhoo 4 and A. Poostforush 1 1 Department of Physics
and Biruni Observatory, Shiraz University, Shiraz 71454, Iran222Permanent
address
2Research Institute for Astronomy and Astrophysics of Maragha, PO Box
55134-441, Maragha, Iran
3 Department of Physics and Astronomy, University of Waterloo, 200 University
Avenue West, Waterloo, Ontario, N2L3G1, Canada
4 Department of Physics, University of Birjand, Birjand, Iran

###### Abstract

In this paper, we investigated the thermodynamic properties of strange quark
matter using Nambu-Jona-Lasinio (NJL) model at finite temperatures where we
considered the dynamical mass as the effective interaction between quarks. By
considering the pressure of strange quark matter (SQM) at finite temperatures,
we showed that the equation of state of this system gets stiffer by increasing
temperature. In addition, we investigated the energy conditions and stability
of the equation of state and showed that the equation of state of SQM satisfy
the conditions of stability. Finally, we computed the structure properties of
hot strange quark stars (SQS) including the gravitational mass, radius,
Schwarzschild radius, average density, compactness and gravitational redshift.
Our calculations showed that in this model, the maximum mass and radius of SQS
increase by increasing temperature. Furthermore it was shown that the average
density of SQS is greater than the normal nuclear density, and it is an
increasing function of temperature. We also discussed the temperature
dependence of the maximum gravitational mass calculated from different
methods.
Keywords: Strange quark matter - strange quark star - NJL model - dynamical
mass - finite temperature

##  I Introduction

By the mid-1970s, physicists realized that hadrons are made up of new
particles later called quarks with a model first proposed by Gell-Mann and
Zwieg Gell-Mann1964 ; Zweig1964 . Baryons at high enough densities (≈$10^{15}$
$gr/cm^{3}$) overlap and dissolve to their components, quarks.

The concept of strange quark matter (SQM), dates back to the works of Jaffe
Jaffe1977 , Chin and Kerman Chin1979 . SQM contains the light quarks (up, down
and strange). In 1984, Witten Witten1984  proposed that SQM might be
absolutely stable, and might be the true ground state of baryonic matter.

Strange quark stars (SQS) are compact objects interesting for astrophysicists
and physicists as the SQS is a great laboratory to study the properties of SQM
due to the density of about $10^{15}$ $g/cm^{3}$. The composition of SQS was
first proposed by Itoh Itoh1970  with the formulation of quantum
chromodynamics (QCD). In 1971, Bodmer Bodmer1971  discussed the possibility of
forming a quark star after the collapse of a massive star, later the concept
of SQS was also mentioned by Witten Witten1984 .

The collapse of a massive star could lead to the formation of a pure SQS by
type $IIa$ supernova (SNII) Sato1987 ; Suzuki1987 ; Hatsuda1987 . Also, a
hybrid star which is a neutron star with a core consisting SQM, can be formed
after neutron star, if the density of the core is high enough. The recent
Chandra observations indicate that objects _RXJ185635-3754_ and _3C58_ may be
SQSs Prakash2003 , as well as candidate for SQS is the object
_SWIFTJ1749.4-2807_ Yu2010 . Actually, a SQS or a hybrid star is denser than a
neutron star. In other words, the mass of SQS is near that of a neutron star
but with a smaller radius.

There are two main frameworks usually used to investigate the thermodynamic
properties of SQM, Nambu-Jona-Lasinio (NJL) model Carroll2009  and MIT bag
model Chodos1974 ; Alford2005 , where theoretical foundations of both is QCD
Freedman1977 . In recent years, we have investigated the thermodynamic
properties of SQM and structure of SQS under different conditions using these
frameworks. We have computed the structural properties of SQS at zero and
finite temperatures, as well as the structure of a magnetized SQS using the
MIT bag model with the fixed and density dependent bag constants at zero and
finite temperatures in the presence and absence of magnetic fields
Bordbar2011vol11 ; Bordbar2011vol54 ; Bordbar2012vol12 ; Bordbar2013 ;
Bordbar2014 . We have also computed the maximum gravitational mass and other
structural properties of a neutron star with a quark core at zero Bordbar2006
and finite temperatures Yazdizadeh2011 .

In our previous work, we have studied the effect of dynamical quark mass in
the calculation of SQS structure using MIT bag model and NJL model at zero
temperature Bordbar2012vol12540 . In current paper, we extend NJL model for
finite temperatures to survey the thermodynamic properties of a hot SQS.
Furthermore, we show that the equation of state of SQM calculated according to
NJL model satisfies the stability and energy conditions. We investigate the
structure properties of SQS by calculating the structure parameters (mass,
radius, Schowarzschild radius, average density, compactness and gravitational
redshift) in the last section.

##  II Calculation of energy and equation of state of hot SQM using NJL model

###  II.1 Nambu-Jona-Lasinio (NJL) model at finite temperatures

The NJL model is named after Nambu, Jona and Lasinio, who for the first time
offered a theory about the dynamical model of elementary particles based on
analogy with super conductivity in $1961$ Nambu1961 . The NJL model is an
effective lagrangian of relativistic fermions interacting through local
fermion-fermion coupling. This model is a suitable approximation of QCD in the
low energy and long wavelength limits Klevansky1992 ; Vogl1991 , appropriate
for the bound states of many-body systems Buballa2005  and EOS of compact
stars Schertler1999 . At high temperatures and densities, interaction leads to
spontaneously breaking of chiral symmetry. In NJL model, symmetry breaking is
characterized by quarks dynamical mass Peng1999 ; Shao2011 .

Dynamical mass of quarks is calculated via,

$M_{i}=m_{0}^{i}-4G<\overline{q_{i}}q_{i}>+2K<\overline{q_{j}}q_{j}><\overline{q_{k}}q_{k}>,$
(1)

where $M_{i}$ is dynamical mass of quark $i$, $m^{i}_{0}$ is mass of free
quark $i$, $G$ and $K$ are the coupling constants, and
$<\overline{q_{i}}q_{i}>$ shows the condensation of quark-antiquark which is
calculated as follows Pennington2005 ; Buballa2005 ,

$<\overline{q_{i}}q_{i}>=-\frac{3}{\pi^{2}(\hbar
c)^{3}}\int_{p_{f}^{i}}^{\lambda}\frac{M_{i}}{\sqrt{M_{i}^{2}c^{2}+p^{2}c^{2}}}f_{i}(p)p^{2}dp.$
(2)

In the above equation, $\lambda$ in the upper limit of integral is the cut-off
value, $p$ is momentum of quark, $p_{f}^{i}$ is the Fermi momentum of each
quark and

$f_{i}(p)=\frac{1}{e^{\beta{(\epsilon_{i}(p)-\mu_{i})}}+1},$ (3)

where $\mu_{i}$ and $\epsilon_{i}$ are the chemical potential and single
particle energy of quark $i$, respectively, and $\beta$ $=$ $1/k_{B}T$
($k_{B}$ is the Boltzmann constant). We calculate the chemical potential,
$\mu_{i}$ by solving the Fermi-Dirac equation numerically. The NJL model is a
renormalizable model, so we should choose a method to find the physical
values. In the present paper, we use an ultra-violent cut-off that indicates
restoring of chiral symmetry breaking, $\lambda=602.3$ $MeV$ Raha2000 ;
Buballa2005 .

###  II.2 Energy and EOS of hot SQM

The total energy density of SQM is defined as the sum of the kinetic energy of
free quarks, $\varepsilon_{i}$, and the potential energy of our system,
$B_{eff}$, which is called the effective bag constant,

$\varepsilon_{tot}=\sum_{i=u,d,s}\varepsilon_{i}+B_{eff},$ (4)

where the kinetic energy of quark $i$ ($\varepsilon_{i}$) is calculated using
the following constraint,

$\varepsilon_{i}=-\frac{3}{\pi^{2}(\hbar
c)^{3}}\int_{0}^{\lambda}\sqrt{M_{i}^{2}c^{2}+p^{2}c^{2}}f_{i}(p)p^{2}dp.$ (5)

The effective bag constant is calculated by the following relation,

$B_{eff}=B_{0}+B_{tot},$ (6)

where,

$B_{tot}=\sum_{i}B_{i}+\frac{1}{(\hbar
c)^{3}}4K<\overline{u}u><\overline{d}d><\overline{s}s>,$ (7)

and,

$B_{i}=\frac{3}{\pi^{2}(\hbar
c)^{3}}\int_{0}^{\lambda}\left[\sqrt{M_{i}^{2}c^{2}+p^{2}c^{2}}-\sqrt{m_{i}^{2}c^{2}+p^{2}c^{2}}-2G<\overline{q_{i}}q_{i}>^{2}\right]p^{2}dp.$
(8)

We need the Helmholtz free energy to calculate the equation of state (EOS) of
the system,

${\cal F}_{tot}=\varepsilon_{tot}-TS_{tot},$ (9)

where $S_{tot}$ is the total entropy of system,

$S_{tot}=\sum_{i=u,d,s}s_{i},$ (10)

and $s_{i}$ is entropy of quark $i$,

$s_{i}=-\frac{3}{\pi^{2}(\hbar
c)^{3}}\int_{0}^{\lambda}\left\\{f_{i}(p)\ln{f_{i}(p)}+[1-f_{i}(p)]\ln{[1-f_{i}(p)]}\right\\}p^{2}dp.$
(11)

To calculate EOS of our system, we use the following relation,

$P(n,T)=\sum_{i}n_{i}\frac{d{\cal F}_{tot}}{dn_{i}}-{\cal F}_{tot},$ (12)

where $n_{i}$ is the number density of quark $i$.

###  II.3 Results of thermodynamic properties of hot SQM

Figure 1: Dynamical mass of $up$ quark versus baryon number density at
different temperature. Figure 2: Dynamical mass of $down$ quark versus baryon
number density at different temperature. Figure 3: Dynamical mass of $strange$
quarks versus baryon number density at different temperature.

In Figs.1, 2 and 3, we have presented the dynamical mass of the up, down and
strange quarks versus baryonic number density, respectively. We compare our
results at different temperatures. As the plots show, the dynamical mass of
each quark tend to the inertial mass ($m_{s}=140.7$ $MeV$ and
$m_{u}=m_{d}=5.5$ $MeV$) by increasing the baryonic density. Also, we can see
for $u$ and $d$ quarks, the dynamical mass increases by increasing the
temperature. These results hold for strange quarks as well, except for density
of $0.5-1\ fm^{-3}$. Our results are also consistent with the previous ones
Bordbar2012vol12540  and the results of Ruster et al. Ruster2006 .

We have shown the total free energy per volume of hot SQM as a function of the
baryonic density in Fig. 4. By increasing the baryonic number density, the
free energy increases. Also it is seen that the free energy decreases with
increasing temperature.

Figure 4: The total free energy per volume of SQM as a function of the
baryonic density at different temperatures.

The pressure of hot SQM at different temperatures has been plotted in Fig. 5.

Figure 5: The pressure of SQM at different temperatures.

This figure shows that the pressure of SQM increases by increasing the
density. We can also see that the pressure increases by increasing the
temperature. These results indicate that the equation of state of SQM becomes
stiffer by increasing the temperature. In other words, the compressibility of
the degenerate gas decreases by increasing temperature, therefore, the EOS
becomes stiffer.

Figure 6: The pressure of SQM versus mass density at different temperatures.

In Fig. 6 we have plotted the pressure of SQM versus mass density at different
temperatures. Our results show that the pressure increases by increasing mass
density. Also, it is shown that the central pressure increases as a function
of temperature.

Here, we show that in the considered version of the equation of state, the
Bodmer-Witten hypothesis holds true. According to this hypothesis, the energy
per particle of SQM should be lower than that of ${}^{56}Fe$ which is
$930.4$$MeV$, so SQM is more stable than the nuclear matter Bodmer1971 ;
Witten1984 . To investigate this condition, we have investigated the energy
per particle behavior at different temperatures ($T$). We have found that the
minimum point of energy per particle versus baryon density which is
corresponding to the zero pressure is equal to $408.77$$MeV$ at $T=30$$MeV$
and is equal to $928.55$$MeV$ at $T=90$$MeV$ which ensures the stability of
SQM. We also study energy and stability conditions in next parts.

###  II.4 Energy conditions

There are four different energy conditions that we study in this work;

a) Null energy condition (NEC) $\longrightarrow$ $P_{c}$ \+
$\rho_{c}c^{2}$$\geq$0,
b) Weak energy condition (WEC) $\longrightarrow$ $P_{c}$ \+
$\rho_{c}c^{2}$$\geq$0 and $\rho_{c}$$\geq$ 0,
c) Strong energy condition (SEC) $\longrightarrow$ $P_{c}$ \+
$\rho_{c}c^{2}$$\geq$ 0 and $3P_{c}$ \+ $\rho_{c}c^{2}$$\geq$ 0,
d) Dominate energy condition (DEC) $\longrightarrow$
$\rho_{c}c^{2}$$\geq$$\mid$ $P_{c}$ $\mid$,

where $\rho_{c}$ and $P_{c}$ are mass density and pressure at the center of
SQS ($r=0$). Results shown in Table 1 at different temperatures correspond to
Fig. 6 and the above four conditions. It is clear that all energy conditions
are satisfied regarding the equation of state we calculated for SQM.

Table 1: Energy conditions of SQS at different temperatures. $T(MeV)$ |  |  $\rho_{c}(10^{15}$$\frac{g}{cm^{3}}$) |  $P_{c}(10^{15}$$\frac{g}{cm^{3}}$) |  | $NEC$ |  | $WEC$ |  | $SEC$ |  | $DEC$
---|---|---|---|---|---|---|---|---|---|---|---
$30$ |  | 7 | 0.969 |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$
$45$ |  | 7 | 1.315 |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$
$70$ |  | 7 | 2.007 |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$
$90$ |  | 7 | 2.630 |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$
$150$ |  | 6.5 | 4.153 |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$ |  | $\checkmark$

###  II.5 Stability of equation of state

To verify the stability of EOS of SQM we use the extreme condition of sound
velocity. The sound velocity is calculated by $v_{s}=\sqrt{dP/d\rho}$. It is
clear that to have a physical model, the sound velocity must satisfy the
condition of $0$ $\leq$$v_{s}^{2}$$\leq$ $c^{2}$. Here we have found that for
all relevant densities and temperatures, the above condition is obeyed by the
velocity of sound. This indicates that the stability of our EOS is confirmed
for all temperatures and densities except densities less than $0.7\times
10^{15}$$gr/cm^{3}$ at temperature of $30$$MeV$. It is clear that the strange
quark matter can be created at high enough temperature and density
Nakazato2008 ; Nakazato2010 ; Nakazato2013 .

##  III Structure properties of strange quark star

The structure of stars is usually determined by their mass and radius,
although there are some other parameters, such as schwarzschild radius,
average density, compactness and gravitational redshift, which we investigate.

###  III.1 Mass and radius of SQS

Since quark stars are relativistic objects we should use the relativistic
equation of hydrostatic equilibrium for these systems,

$\frac{dP}{dr}=-\frac{G\left[\varepsilon(r)+\frac{P(r)}{c^{2}}\right]\left[m(r)+\frac{4\pi
r^{3}P(r)}{c^{2}}\right]}{r^{2}\left[1-\frac{2Gm(r)}{rc^{2}}\right]},$ (13)
$\frac{dm}{dr}=4\pi r^{2}\varepsilon(r).$ (14)

These equations are known as Tolman-Oppenheimer-Volkov equations (TOV)
Oppenheimer1939 . Using the equation of state which was obtained in the
previous section and the boundery conditions ($P(r=0)=P_{c}$, $P(r=R)=0$,
$m(r=0)=0$ and $m(r=R)=M_{max}$) we integrate the TOV equations to compute the
structure of strange quark stars (SQS).

Figure 7: The gravitational mass of SQS versus energy density at different
temperatures.

In Fig. 7, we have plotted the gravitational mass of strange quark star (SQS)
versus energy density at different temperatures. We can see that for all
temperatures, the gravitational mass increases rapidly by increasing the
energy density, and finally reaches a limiting value (maximum gravitational
mass). The maximum gravitational mass for different temperatures has been
given in Table 2. Our results show that this maximum mass increases by
increasing the temperature. We have shown the gravitational mass of SQS as a
function of the radius (M-R relation) at different temperatures in Fig. 8.

Figure 8: The gravitational mass of SQS as a function of the radius at
different temperatures.

This figure shows that by increasing the gravitational mass, till the maximum
mass is reached, the radius increases. We can see that the increasing rate of
gravitational mass versus radius increases by increasing temperature. The
radius of SQS corresponding to the maximum mass has been given in Table 2
indicating higher radius for higher temperatures. Here, it should be noted
that as it was seen from Fig. 5, by increasing the temperature, the equation
of state of SQS becomes stiffer. Now, we can conclude that in the finite
temperature NJL model of SQS, the stiffer equation of state leads to the
higher maximum gravitational mass for this compact object (Table 2). This
behavior has been also reported by Chu et al. Chu2017 .

Table 2: Structure properties of SQS at different temperatures. $T(MeV)$ | $M_{max}(M_{\odot})$ | $R(km)$ |  | $\overline{\rm\rho}(10^{15}\frac{g}{cm^{3}})$ |  | $\sigma$ |  | $Z_{s}$
---|---|---|---|---|---|---|---|---
$30$ | 0.650 | 5.305 |  | 1.86 |  | 0.348 |  | 0.239
$45$ | 0.742 | 5.498 |  | 2.29 |  | 0.411 |  | 0.303
$70$ | 0.85 | 5.962 |  | 2.12 |  | 0.436 |  | 0.331
$90$ | 1.002 | 6.339 |  | 1.85 |  | 0.448 |  | 0.346

###  III.2 Average density

We can calculate the average density of the star using maximum mass ($M$) and
radius ($R$) by,

$\overline{\rm\rho}=\frac{3M}{4\pi R^{3}}.$ (15)

The results of this calculation are shown in Table 2. The minimum average
density regarding Table 2, $\overline{\rm\rho}=1.85\times 10^{15}g/cm^{3}$, is
related to temperature of $90MeV$ which is larger than the normal nuclear
matter density, $\rho_{0}=2.7\times 10^{14}g/cm^{3}$. Furthermore, the central
density of SQS regarding Table 1 is about $7\times 10^{15}g/cm^{3}$, which is
larger than average density of SQS at all temperatures.

###  III.3 Compactness

The compactness is a parameter to show the strength of gravity. It is
calculated using the ratio of Schwarzschild radius to radius of star
($\sigma=R_{sch}/R$ where $R_{sch}=\frac{2GM}{c^{2}}$). As it is shown in
Table 2, $\sigma$ is almost the same at all temperatures for SQS.

###  III.4 Gravitational redshift

The gravitational redshift is calculated as,

$Z_{s}=\frac{1}{\sqrt{1-\frac{2GM}{c^{2}R}}}-1,$ (16)

where $M$ is the maximum mass and $R$ is the radius of SQS. We have plotted
the gravitational redshift of SQS versus the gravitational mass at different
temperatures in Fig. 9.

Figure 9: The gravitational redshift of SQS as a function of the gravitational
mass at different temperatures.

Obviously, it can be seen that for all temperatures, the gravitational
redshift increases by increasing the gravitational mass to the value of the
maximum limit. Also, it is clear that the gravitational redshift increases by
increasing temperature. The results of gravitational redshift of SQS
corresponding to maximum mass and temperature have been shown in the last
column of Table 2. It can be seen that the gravitational redshift increases as
a function of maximum mass. The maximum gravitational redshift is calculated,
$z_{s}$=$0.346$, at temperature of $T=90$ $MeV$, that is about $59.34$$\%$
less than critical value of gravitational redshift ($Z_{s}^{CL}$=$0.8509$)
Haensel2007 . Furthermore, the gravitational redshift at temperature of
$T$=$90$$MeV$ ($Z_{s}$=$0.346$) is about $0.4$$\%$ less than the observational
result that is reported for quark star candidate $RXJ185635-3754$
($Z_{s}$=$0.35$ $\pm$ $0.15$) Prakash2003 .

###  III.5 The mass of SQS in term of Planck mass

In this section, we show that the mass of SQS can be expressed in term of the
fundamental value of Planck mass, then we derive the relevant relation. The
repulsive nuclear force and the degeneracy pressure of fermions both are
against gravity to avoid the collapse of compact stars. As we have mentioned,
by phase transition of nucleons, the density of SQS is near and above the
normal nuclear matter density. Therefore, using these facts, in the maximum
value, we can consider the average density of SQS equal to the nuclear
density, where the nuclear density is approximately defined as follow,

$\rho_{nuc}\simeq\frac{3m_{p}}{4\pi\lambda_{\pi}^{3}},$ (17)

where $m_{p}$ is the proton mass and $\lambda_{\pi}=\hbar/m_{\pi}c$ is the
Compton wavelength of pion. From previous sections, we use $R_{sch}$ and
$\overline{\rm\rho}$ to derive the following equation,

$M\simeq(\frac{\hbar
c}{G})^{3/2}\frac{1}{m_{p}^{2}}(\frac{\eta_{\pi}}{2\eta_{p}})^{3/2}\simeq
M_{Ch}(\frac{\eta_{\pi}}{2\eta_{p}})^{3/2}\simeq
m_{pl}\eta_{p}^{2}(\frac{\eta_{\pi}}{2\eta_{p}})^{3/2},$ (18)

where $M_{Ch}$ is the Chandrasekhar mass ($\simeq(\frac{\hbar
c}{G})^{3/2}\frac{1}{m_{p}^{2}}$), $m_{pl}$ is Planck mass,
$\eta_{p}=m_{pl}/m_{p}$ and $\eta_{\pi}=m_{pl}/m_{\pi}$ Burrows2014 .

##  IV  The temperature dependence of graviational maximum mass of SQS

In this section, we want to look at the behavior of the maximum gravitational
mass of SQS which is calculated from different methods at finite temperature.

We have calculated the thermodynamic properties and structure of SQS at finite
temperature using MIT bag model with the fixed bag constant and density-
dependent bag constant Bordbar2011vol54 . It has been shown that the EOS of
the system in both cases (fixed bag constant and density-dependent bag
constant) becomes stiffer by increasing temperature. Then, we have shown that
the maximum gravitational mass and the corresponding radius decrease as a
function of temperature in both mentioned cases. For $B=90$$MeV$, the maximum
gravitational mass and the corresponding radius have been calculated to be
$1.228$ $M_{\odot}$ and $7.073$ $km$ at $T=30$ $MeV$ and $1.04$ $M_{\odot}$
and $6.14$ $km$ at $T=80$ $MeV$. Also for density dependent bag constant the
maximum gravitational mass has been changed from $1.34$ to $1.12$ $M_{\odot}$
where temperature changed from $T=30$ $MeV$ to $T=80$ $MeV$. In the same way
the radius of SQS decreases from $7.44$ to $6.57$ $km$.

We have also investigated the structure of spin-polarized strange quark star
at finite temperature using MIT bag model with $B=90$$MeV$ Bordbar2011vol11 ,
and with a density-dependent bag constant Bordbar2012vol12 . The EOS and the
maximum gravitational mass and radius in Refs. Bordbar2013 ; Bordbar2014  were
similar to the previous work Bordbar2011vol54 . The maximum gravitational mass
decrease from $1.171$ $M_{\odot}$ at $T=30$ $MeV$ to $1.16$ $M_{\odot}$ at
$T=70$ $MeV$ using $B=90$ $MeV$, and the radius decreased from $7.27$ to
$7.21$ $km$. By considering a density dependent bag constant the maximum
gravitational mass and the corresponding radius decreased from
$1.15$$M_{\odot}$ and $7.1$ $km$ at $T=30$ $MeV$ to $0.77$ $M_{\odot}$ and
$6.89$ $km$ at $T=70$ $MeV$.

The structure of SQS has been calculated by Alaverdyan and Hajyan
Alaverdyan2014 . They have considered the ultrarelativistic quarks in SQS and
have calculated the EOS of the system using MIT bag model. As they have
reported, the EOS becomes stiffer as the temperature increases, where the
radius of SQS versus temperature has been plotted. It can be seen from this
figure that the radius and the corresponding gravitational mass increase from
$7.23$ $km$ and $0.49$ $M_{\odot}$ to $8.27$ $km$ and $0.77$ $M_{\odot}$ when
temperature increases in a range from zero to $80$ $MeV$.

Compact strange stars with a medium dependence on gluons at finite temperature
have been studied by Bagchi et al. Bagchi2006 . Properties have been
calculated using large color approximation with built-in chiral symmetry
restoration in that paper. Their calculations have shown that the stiffer EOS
has been achieved at higher temperatures. Similarity, the maximum
gravitational mass and the correspondig radius are larger at the lower
temperatures.

As we have shown in the section III.1 of the present paper, using NJL model
creates a different behavior in gravitational mass as a function of
temperature. We can see from Table 2, the maximum gravitational mass and
radius increase by increasing temperature, although the EOS of the system
becomes stiffer by increasing temperature.

A similar behavior with our current work also has been reported in Ref.
Chu2017 . They have used NJL model as in our current paper. They have plotted
the equation of state at three different temperatures. It has been shown that
EOS becomes stiffer by increasing temperature as we have shown in Fig. 5.
Furthermore, the gravitational mass as a function of temperature has been
plotted and it has been reported that when temperature rises to $50$, $80$,
and $100$ $MeV$, the maximum mass of quark stars will reach $2.13$
$M_{\odot}$, $2.46$ $M_{\odot}$, and $2.71$ $M_{\odot}$, respectively.

In another work Dexheimer2013  the authors have compared the properties of the
proto-quark star from different methods (Quark-mass density-dependent (QMDD)
model, MIT bag model, and NJL model) at finite temperature. They have reported
the maximum gravitational mass which is calculated from QMDD and MIT bag
models at different temperatures. Their results show that investigation of the
proto-quark star gives different behavior for mass-radius results as the
function of temperature by MIT bag model. Where they have considered the same
conditions as in Ref. Bordbar2011vol54 , different behavior has been achieved
for mass and radius in the temperature range. The gravitational mass increases
from $1.62$ to $1.65$ $M_{\odot}$ and the corresponding radius increases from
$9.01$ to $9.15$ $km$. So, when the conditions are such as in Refs.
Bordbar2013 ; Bordbar2014 , although the gravitational mass decreases from
$2.02$ to $1.93$ $M_{\odot}$, the radius increases from $9.04$ to $9.08$ $km$
by increasing temperature. Using QMDD, they have considered two versions: $1)$
where the pressure at the density corresponding to the minimum of the free
energy per baryon could be non-zero, depending on the matter studied (SM or
2QM) is noted as version 1 (QMDDv1) $2)$ presenting a remedy to the
thermodynamical inconsistency, in such a way that the minimum of the energy
per baryon corresponds to the point of zero pressure, is noted as version 2
(QMDDv2). In both versions, they have considered different masses for strange
quarks ($150$ and $100$ $MeV/c^{2}$ ). Using QMDDv1, it has been shown that
the maximum gravitational mass and the corresponding radius increase by
increasing temperature for both strange quark masses (from $2.28$ to $2.33$
$M_{\odot}$ and $12.05$ to $12.19$ $km$ for the first strange quark mass and
from $2.26$ to $2.29$ $M_{\odot}$ and $11.75$ to $11.76$ $km$ for the second
strange quark mass). Using QMDDv2 when they have used strange quark mass equal
to $150$ $MeV/c^{2}$ the maximum gravitational mass and corresponding radius
increases by increasing temperature (from $1.60$ to $1.62$ $M_{\odot}$ and
$8.42$ to $8.46$ $km$), but when they have considered the strange quark mass
equal to $100$ $MeV/c^{2}$ there is an inverse behavior for the maximum
gravitational mass and radius in the temperature range (the gravitational mass
decrease from $1.59$ to $1.58$ $M_{\odot}$ and the radius decrease from $8.22$
to $8.16$ $km$).

##  V Summary and conclusions

In this paper, we have calculated the thermodynamic properties of the strange
quark matter (SQM) at finite temperatures using NJL model, and we have
investigated the structure of strange quark stars (SQS). We have calculated
free energy and equation of state (EOS) of SQM by considering the dynamical
mass. We have shown that free energy increases the corresponding baryonic
density. In addition, free energy decreases by increasing temperature at a
specified density. Also, our results indicate that pressure increases
proportional to the density and EOS of SQM becomes stiffer as a function of
temperature. Furthermore, we have investigated the energy conditions and
stability of EOS. We have shown that EOS of our system satisfies both energy
conditions and stability.

Later, we studied the structure of SQS using the general relativistic TOV
equations and boundary conditions. We calculated the maximum gravitational
mass and the corresponding radius of SQS at different temperatures. Follow up
on the structure of the star, we have calculated other parameters such as the
Schwarzschild Radius, average density of SQS, compactness and gravitational
redshift. We have shown that the gravitational mass and radius of SQS rapidly
increases by increasing temperature and we have compared the behavior of
temperature dependent maximum gravitational mass for different methods. We
have shown that the average density of SQS is more than the normal nuclear
matter density. In addition, our calculations show that the compactness of SQS
is almost the same at all temperatures. We have also investigated the
gravitational redshift ($Z_{s}$) of SQS at different temperatures and found
that $Z_{s}$ increases by increasing temperature. Comparison with the
observational results clarify that the gravitational redshift of SQS at
temperature of $90$$MeV$ is just $\%$$3$ less than the gravitational redshift
of $RXJ185635-3754$. Finally we have derived the relation between mass of SQS
and Plank mass.

## Acknowledgements

We wish to thank Shiraz University Research Council. This work has been
supported by Research Institute for Astronomy and Astrophysics of Maragha. G.
H. Bordbar wishes to thank A. Broderick (University of Waterloo) for his
useful comments and discussions during this work. G. H. Bordbar also wishes to
thank Physics Department of University of Waterloo for the great hospitality
during his sabbatical.

## References

  * (1) M. Gell-Mann, _Phys. Lett._ 8, 214 (1964).
  * (2) G. Zweig, _Cern-Reports, TH-401, TH-412_ (1964).
  * (3) R. L. Jaffe, _Phys. Rev. Lett._ 38, 195, 617E (1977).
  * (4) S. A. Chin and A.K. Kerman, _Phys. Rev. Lett._ 43, 1292 (1979).
  * (5) E. Witten, _Phys. Rev._ D30, 272 (1984).
  * (6) N. Itoh, _Prog. Theor. Phys._ 44, 291 (1970).
  * (7) A. R. Bodmer, _Phys. Rev._ D4, 1601 (1971).
  * (8) K. Sato and H. Suzuki, _Phys. Rev. Lett_ 58, 2722 (1987).
  * (9) H. Suzuki and K. Sato, _Perprint. UTAP_ , 53/87 (1987).
  * (10) T. Hatsuda, _Mod. Phys. Lett._ A2, 805 (1987).
  * (11) M. Prakash, J. M. Lattimer, A. W. Steiner and D. Page, _Nucl. Phys._ A715, 835 (2003).
  * (12) H. W. Yu and R. X. Xu, _Res. Astron. Astrophys._ 11, 471 (2010).
  * (13) J. D. Carroll, D. B. Leinweber, A. W. Thomas, and A. G. Williams, _Phys. Rev._ C79, 045810 (2009).
  * (14) A. Chodos, R. L. Jaffe, K. Johnson, C. B. Thorn and V. F. Weisskopf, _Phys. Rev._ D9, 3471 (1974).
  * (15) M. Alford, M. Braby, M. Paris and S. Reddy, _Astrophys. J._ 626, 969 (2005).
  * (16) B. Freedman and L. Mclerran, _Phys. Rev._ D16, 1130 (1977).
  * (17) G. H. Bordbar and A. Peivand, _Res. Astron. Astrophys._ 11, 851 (2011).
  * (18) G. H. Bordbar, A. Poostforush and A. Zamani, _Astrophysics_ , 54, 277 (2011).
  * (19) G.H. Bordbar, H. Bahri and F. Kayanikhoo, _Res. Astron. Astrophys._ 12, 1280 (2012).
  * (20) G. H. Bordbar, F. Kayanikhoo and H. Bahri, _Iranian J. Sci. Tech._ A37, 165 (2013).
  * (21) G. H. Bordbar and Z. Alizadeh, _Astrophysics_ 57, 130 (2014).
  * (22) G. H. Bordbar, M. Bigdeli and T. Yazdizadeh, _Int. J. Mod. Phys._ A21, 5991 (2006).
  * (23) T. Yazdizadeh and G. H. Bordbar, _Res. Astron. Astrophys._ 11, 471 (2011).
  * (24) G. H. Bordbar and B. Ziaei, _Res. Astron. Astrophys._ 12, 540 (2012).
  * (25) Y. Nambu and G. Jona-Lasinio, _Phys. Rev._ 122, 345 (1961).
  * (26) S. P. Klevanski, _Rev. Mod. Phys._ 64, 3 (1992).
  * (27) U. Vogl and W. Weise, _Prog. Part. Nucl. Phys._ 27, 195 (1991).
  * (28) M. Buballa, _Phys. Rep._ 407, 205 (2005).
  * (29) K. Schertler, S. Leupold, and J. Schaffner-Bielich, _Phys. Rev._ C60, 025801 (1999).
  * (30) G. X. Peng, H. C. Chiang, J. J. Yang, L. Li, and B. Liu, _Phys. Rev._ C61, 015201 (1999).
  * (31) G. Y. Shao, M. Di Toro, B. Liu, M. Colonna, V. Greco, Y. X. Liu, and S. Plumari, _Phys. Rev._ D83, 094033 (2011).
  * (32) M. R. Pennington, _J. Phys. Conf. Ser._ 18, 1 (2005).
  * (33) S. Raha, _AIP Conference Proceedings_ 508, 226 (2000).
  * (34) S. B. R¨uster, V. Werth, M. Buballa, I. A. Shovkovy and D. H. Rischke, _Phys. Rev._ D73, 034025 (2006).
  * (35) K. Nakazato, K. Sumiyosh and S. Yamada, _Phys. Rev._ D77, 103006, (2008).
  * (36) K. Nakazato, K. Sumiyosh and S. Yamada, _Astrophys. J._ 721, 1284, (2010).
  * (37) K. Nakazato, K. Sumiyosh and S. Yamada, _Astron.Astrophys_ A50, 558, (2013).
  * (38) J. R. Oppenheimer and G. M. Volkoff, _Phys. Rev._ 55, 374 (1939).
  * (39) P. Chu, X. Li, B. Wang, Y. Dong, Y. Jia, S. Wang, H. Ma, _Eur. Phys. J._ C77, 512 (2017).
  * (40) P. Haensel, A. Y. Potekhin, and D.G. Yakovlev, Neutron stars 1: Equation of state and structure, _Springer_ (2007).
  * (41) A. Burrows and J. P. Ostriker, _Proc. Nat. Acad. Sci._ 111, 2409 (2014).
  * (42) A. G. Alaverdyan, G. S. Hajyan, _J. Phys.: Conf. Ser._ 496, 012005 (2014).
  * (43) M. Bagchi, S. Ray, M. Dey and J. Dey, _Astron. Astrophys._ 450, 431 (2006).
  * (44) V. Dexheimer, J. R. Torres and D. P. Menezes, _Eur. Phys. J._ C73, 2569 (2013).


# A Refined Proximal Algorithm for Nonconvex Multiobjective Optimization in
Hilbert Spaces

G. C. Bento∗ ,  J. X. Cruz Neto  ,  J. O. Lopes  ,  B. S. Mordukhovich  and
P. R. Silva Filho  Institute of Mathematics and Statistics, IME, Federal
University of Goiás, Goiânia, Goiás, Brazil<EMAIL_ADDRESS>Departamento de
Matemática, CCN, Federal University of Piauí
Teresina, Piauí, Brazil<EMAIL_ADDRESS>Departamento de Matemática, CCN,
Federal University of Piauí
Teresina, Piauí, Brazil<EMAIL_ADDRESS>Department of Mathematics,
Wayne State University, Detroit, Michigan, United States<EMAIL_ADDRESS>Departamento de Matemática, CCN, Federal University of Piauí
Teresina, Piauí, Brazil<EMAIL_ADDRESS>
###### Abstract.

This paper is devoted to general nonconvex problems of multiobjective
optimization in Hilbert spaces. Based on Mordukhovich’s limiting subgradients,
we define a new notion of Pareto critical points for such problems, establish
necessary optimality conditions for them, and then employ these conditions to
develop a refined version of the vectorial proximal point algorithm with
providing its detailed convergence analysis. The obtained results largely
extend those initiated by Bonnel, Iusem and Svaiter [4] for convex vector
optimization problems and by Bento et al. [3] for nonconvex finite-dimensional
problems in terms of Clarke’s generalized gradients.

###### Key words and phrases:

Multiobjective programming, proximal point algorithm, Locally Lipschitz
functions, weakly Pareto point

##  1\. Introduction

Multiobjective optimization involves optimizing two or more real-valued
objective functions simultaneously. Typically, no single point can minimize
all the given objective functions simultaneously, which leads us to the
concept of Pareto optimality replacing the traditional notion of optimality.
Restricting our discussion to the proximal point algorithm (PPA), the focus of
this paper is to improve the (exact) method introduced and analyzed by Bonnel,
Iusem sand Svaiter [4]. In the mentioned paper, the authors employ the method,
in the fully convex setting, to find a weak Pareto minimizer of a vector
function $F:X\to Y$ from a real Hilbert space $X$ to a real Banach space $Y$
containing a closed, convex, and pointed cone $\mathcal{C}$ with nonempty
interior. The method constructs as $(k+1)$th iteration a weak Pareto solution
of $F_{k}(x):=F(x)+\lambda_{k}||x-x^{k}||^{2}\varepsilon^{k}$ constrained to
the set $\Omega_{k}:=\\{x\in X\;:\;F(x)\preceq_{\mathcal{C}}F(x^{k})\\}$. Here
$\\{\lambda_{k}\\}$ represents a bounded sequence of positive scalars, and
$\varepsilon^{k}$ denotes an exogenously selected vector from the interior of
$\mathcal{C}$ with $||\varepsilon^{k}||=1$ for each $k=0,1,\ldots$. The
strategy used in the convergence results is based on the first-order
optimality condition for the following scalar problem:

(1.1) $\min_{x\in\Omega_{k}}\eta_{k}(x),$

where $\eta_{k}(x):=\langle
F(x),z^{k}\rangle+\frac{\lambda_{k}}{2}\langle\varepsilon^{k},z^{k}\rangle||x-x^{k}||^{2}$
and $\\{z^{k}\\}\in\mathcal{C}^{+}\subset Y^{*}$, $||z^{k}||=1$, is an
exogenous sequence, $\mathcal{C}^{+}$ denotes the positive polar cone and
$Y^{*}$ is the topological dual space of $Y$ with
$\langle\cdot,\cdot\rangle:Y\times Y^{*}\to\mathbb{R}$ being the duality
pairing. Therefore, solving (1.1) implies that

(1.2)
$0\in\partial\psi_{k}(x^{k+1})+\lambda_{k}\langle\varepsilon^{k},z^{k}\rangle(x^{k+1}-x^{k}),$

where $\psi_{k}(x):=\langle F(x),z^{k}\rangle+\delta_{\Omega_{k}}(x)$, where
$\partial\psi_{k}$ denotes the subdifferential of $\psi_{k}(\cdot)$ in the
sense of convex analysis, and where $\delta_{\Omega_{k}}(\cdot)$ is the
indicator function of $\Omega_{k}$. Bonnel, Iudsem and Svaiter [4] prove that
any sequence generated by the PPA converges (in the weak topology of $X$) to a
weakly efficient minimizer of $F$, provided that a certain completeness
condition of $Y$ is satisfied. Note that (1.2) can be viewed as

$\alpha_{k}(x^{k}-x^{k+1})\in\partial(\langle
F(\cdot),z^{k}\rangle)(x^{k+1})+N(x^{k+1};\Omega_{k}),$

where $\alpha_{k}=\lambda_{k}\langle\varepsilon^{k},z^{k}\rangle$ and
$N(x^{k+1};\Omega_{k})$ stands for the normal cone to $\Omega_{k}$ at
$x^{k+1}\in\Omega_{k}$ in the classical sense of convex analysis. In this
approach, the set $\Omega_{k}$ ensures that the algorithm is functioning as a
descent process.

A motivation, in a dynamic context, to consider the constrained set
$\Omega_{k}$ is provided by Bento et al. [2]. They observe that the set
$\Omega_{k}$ delineates a vectorial improvement process, where a vector
minimizing solution $x^{k}$ of the current proximal problem moves to a next
one in a manner that it improves the current solution. This aspect is crucial
for justifying the process at a behavioral level, where a risk-averse agent
agrees to make changes only if they lead us to improvements across all aspects
(all components of the vector). In the case of a group of agents, as
considered here, this constraint set becomes even more significant. It ensures
that the payoff of each agent in the group does not diminish. For variations
of the algorithm initially introduced in [4] for convex vector problems see,
for example, Ceng and Yao [7], Choung et al. [8], Gregório and Oliveira [11]
and Villacorta and Oliveira [18]. Discussions about the
$\mathbb{R}^{m}_{+}$-quasiconvex case have been explored in the works of Bento
et al.[2] and Apolinario et al.[1]. In these studies, their respective
algorithms compute a point $x^{k+1}$ at the $(k+1)$th iteration satisfying the
equation:

$0\in\partial g(F(x^{k+1}))+\alpha_{k}(x^{k+1}-x^{k})+N(x^{k+1};\Omega_{k}),$

where $g:\mathbb{R}^{m}\to\mathbb{R}$ is a kind of “scalarization” function,
where $\partial g$ denotes some subdifferential of $g$, and where
$\\{\alpha_{k}\\}$ is a sequence of positive numbers. For an alternative
approach that does not necessarily depend on the convexity of $\Omega_{k}$
(the previously mentioned works depend on the convexity of $\Omega_{k}$), we
refer to Bento et al. [3]. In the aforementioned paper, the authors extend the
applicability of the PPA (in finite dimensions) to locally Lipschitz vector-
valued functions. It is worth pointing out that the key aspect of this
alternative approach lies in selecting each iteration of the algorithm as a
weak Pareto solution obtained through an optimality condition directly applied
to the vector subproblem restricted to $\Omega_{k}$; namely, extended
conditions of the Fritz-John type given by Clarke’s generalized gradients,
which are necessary for weak Pareto-optimal solutions to multiobjective
problems.

The aim of this paper is to investigate generalized conditions of the Fritz-
John type in the fully nonconvex setting, characterized by Mordukhovich’s
limiting subgradients, and to broaden the applicability of the PPA (in Hilbert
spaces) to locally Lipschitz vector-valued functions. More specifically, we
aim to ensure that each weak cluster point represents a critical Pareto point
of the problem defined via limiting subgradients. Additionally, under the
supplementary assumptions of pseudoconvexity of the vector-valued function and
completeness of its counter-domain, we demonstrate that thge entire sequence
of iterates generated by the PPA converges (in the weak topology of the space)
to a Pareto critical point of the problem.

The rest of the paper is organized as follows. Section 2 introduces major
notations together with some concepts and results of multiobjective
optimization and generalized differentiation in variational analysis needed
below. In Section 3, we define the proximal point method to solve nonconvex
multiobjective optimization problems and then establish its important
properties. Finally, Section 4 is devoted to the convergence analysis of the
algorithm.

##  2\. Multiobjective optimization and generalized differentiation

In this section, we revisit fundamental definitions and properties in
multiobjective optimization as outlined, for example, in the book by Luc [13].
Furthermore, some notions and results of variational analysis and generalized
differentiation, taken mostly from the books by Mordukhovich [15] and
Rockafellar-Wets [17] and playing a crucial role in the main achievements of
this paper, are also presented here.

We consider the $m$-dimensional Euclidean space $\mathbb{R}^{m}$ with the
partial order denoted by $\preceq$ induced by the Pareto cone
$\mathbb{R}^{m}_{+}$ as follows: $y\preceq z$ (or $z\succeq y$) if and only if
$z-y\in\mathbb{R}^{m}_{+}$. The associated notion $y\prec Z$ means that
$z-y\in\mathbb{R}^{m}_{++}$, where $\mathbb{R}^{m}_{+}$ and
$\mathbb{R}^{m}_{++}$ are defined by

$\mathbb{R}^{m}_{+}:=\left\\{x\in\mathbb{R}^{m}~{}:~{}x_{j}\geq
0,~{}j\in\mathcal{I}\right\\},\quad\mathbb{R}^{m}_{++}:=\left\\{x\in\mathbb{R}^{m}~{}:~{}x_{j}>0,~{}j\in\mathcal{I}\right\\},$

and where $\mathcal{I}:=\\{1,\ldots,m\\}$. Let $\mathcal{H}$ be a real Hilbert
space with the inner product and the norm denoted, respectively, by
$\langle\cdot,\cdot\rangle$ and $\|\cdot\|$. Given a nonempty closed set
$\Omega\subset\mathcal{H}$ and a vector-valued function
$F=(f_{1},\ldots,f_{m}):\mathcal{H}\to\mathbb{R}^{m}$, the problem of finding
a weak Pareto point of $F$ in $\Omega$ consists of determining a point
$x^{*}\in\Omega$ such that there exists no vector $x\in\Omega$ with $F(x)\prec
F(x^{*})$. We denote this problem by

(2.1) $\mbox{min}_{w}\\{F(x)~{}:~{}x\in\Omega\\}.$

By $\mbox{arg}\min_{w}\\{F(x):x\in\Omega\\}$ denote the collection of all weak
Pareto points of the vector function $F$ in the set $\Omega$.

######  Remark 1.

As mentioned in Huang and Yang [12], the vector functions

$F(\cdot)\quad\mbox{and}\quad
e^{F(\cdot)}:=(e^{f_{1}(\cdot)},\ldots,e^{f_{m}(\cdot)})$

have the same sets of weak Pareto points, where $e^{\alpha}$ denotes the
exponential map valued at $\alpha\in\mathbb{R}$. This result can be easily
extended to Pareto critical setting in Definition 1. Therefore, concerning
Pareto critical points, we can assume without loss of generality that $F\succ
0$.

Let $f:\mathcal{H}\to\mathbb{R}$ be a lower semicontinuous function. For a
given sequence $\\{u_{k}\\}\subset\mathcal{H}$, we use the notation
$u^{k}\stackrel{{\scriptstyle f}}{{\to}}\bar{u}$ to signify that
$u^{k}\to\bar{u}$ with $f(u^{k})\to f(\bar{u})$. As mentioned in the
introduction, the objective function is not necessarily differentiable. To
handle such situations, we introduce the following concepts. The
regular/Fréchet subdifferential, of $f$ at $\bar{a}\in\mathcal{H}$ is defined
by

$\partial^{F}f(\overline{a})=\bigg{\\{}v\in\mathcal{H}:\,\liminf_{a\to\overline{a}}\dfrac{f(a)-f(\overline{a})-\langle
v,a-\overline{a}\rangle}{\|a-\overline{a}\|}\geq 0\bigg{\\}}.$

It can be equivalently described as is the set of vectors $v\in\mathcal{H}$
such that

$\displaystyle f(a)\geq f(\bar{a})+\langle
v,a-\overline{a}\rangle+o(\|a-\overline{a}\|),$

where
$\displaystyle\lim_{a\to\overline{a}}\dfrac{o(\|a-\overline{a}\|)}{\|a-\overline{a}\|}=0$.
The limiting/Mordukhovich subdifferential of $f$ at $\bar{a}\in\mathcal{H}$ is
defined by

$\partial^{M}f(\bar{a})=\operatorname*{Lim\,sup}_{a\stackrel{{\scriptstyle
f}}{{\to}}\bar{a}}\partial^{F}f(a)=\left\\{v\in\mathcal{H}:\exists\,a^{k}\stackrel{{\scriptstyle
f}}{{\to}}\bar{a},\,v^{k}\in\partial^{F}f(a^{k}),\,v^{k}\rightharpoonup
v\right\\},$

where “Limsup” signifies the outer limit of set-valued mappings in the norm-
weak topology of ${\mathcal{H}}\times{\mathcal{H}}$.

It is well known that if $f$ is differentiable at $\bar{a}\in\mathcal{H}$,
then its Fréchet subdifferential is a singleton,
$\partial^{F}f(\bar{a})=\\{\nabla f(\bar{x})\\}$. If $f$ is continuously
differentiable around $\bar{a}\in\mathcal{H}$ (or merely strictly
differentiable at this point), then its Mordukhovich subdifferential reduces
to the same singleton.

Let $f:\mathcal{H}\to\mathbb{R}$ be a locally Lipschitz function around
$x\in\mathcal{H}$, and let $d\in\mathcal{H}$. The (Clarke) generalized
directional derivative of $f$ at $x$ in the direction of $d$, denoted by
$f^{\circ}(x,d)$, is defined by

$\displaystyle f^{\circ}(x,d)=\limsup_{{y\to x}\atop{t\downarrow
0}}\dfrac{f(y+td)-f(y)}{t}.$

The Clarke subdifferential/generalized gradient of $f$ at $x$, denoted by
$\partial^{C}f(x)$, is defined via the generalized directional derivative by

$\displaystyle\partial^{C}f(x)=\\{b\in\mathcal{H}:\langle b,d\rangle\leq
f^{\circ}(x,d)\;\mbox{ for all }\;d\in\mathcal{H}\\}.$

Note that, in the Hilbert space setting, we have the following relationship
between the Clarke and Mordukhovich subdifferentials (see [15, Theorem 3.57])

(2.2) $\partial^{C}f(x)={\rm clco}\,\partial^{M}f(x),$

where “clco” stands for the weak closure of the convex hull of the set.

Define further the distance function $d:\mathcal{H}\to\mathbb{R}$ of
$x\in\mathcal{H}$ to $\Omega\subset\mathcal{H}$ by

$d_{\Omega}(x):=\inf\\{||x-y||\;:\;y\in\Omega\\}.$

Recall also that by a lower semicontinuous vector function
$F:\mathcal{H}\to\mathbb{R}^{m}$, we mean that component
$f_{i}:\mathcal{H}\to\mathbb{R}$ of $F$ are lower semicontinuous.

In our study, we will mostly use the limiting subdifferential $\partial^{M}$
whose graph is closed in rather general settings; see below. The regular
subdifferential does not have the closed graph even for simple Lipschitz
functions of one variable; see, e.g., $f(x):=-|x|$, where
$\partial^{F}f(0)=\emptyset$. The usage of Clarke’s subdifferential brings us
to more restrictive results due the convexification in (2.2).

To begin with, we give a new definition of Pareto critical points for
nonconvex problems, which is the limiting subdifferential extension of the
classical notion for convex problems of multiobjective optimization.

######  Definition 1.

We say that $\bar{x}\in\Omega\subset\mathcal{H}$ is a Pareto critical point of
the optimization problem (2.1) if there exist $i\in\mathcal{I}$ and
$\omega_{i}\in\partial^{M}f_{i}(\bar{x})$ such that

$\displaystyle\langle\omega_{i},y-\bar{x}\rangle\geq 0\;\mbox{ for all
}\;y\in\Omega.$

To proceed with the corresponding geometric constructions, the regular/Fréchet
normal cone to a nonempty set $\Omega$ at $x\in\Omega$ is given by

$\displaystyle
N^{F}(x;\Omega):=\left\\{\omega\in\mathcal{H}:\;\limsup_{y\xrightarrow{\Omega}x}\dfrac{\langle\omega,y-x\rangle}{\|y-x\|}\leq
0\right\\},$

where the notation $y\xrightarrow{\Omega}x$ means $y\to x$ with $y\in\Omega$.
The limiting/Mordukhovich normal cone to $\Omega$ at $x\in\Omega$ is defined
by

$\displaystyle N^{M}(x;\Omega)$ $\displaystyle:=\operatorname*{Lim\,sup}_{y\to
x}N^{F}(y;\Omega)$ $\displaystyle=\\{\omega\in\mathcal{H}:\,\exists\
x^{k}\xrightarrow{\Omega}x,\exists\ \omega_{k}\in N^{F}(x_{k};\Omega);\
\omega_{k}\rightharpoonup\omega\\}.$

It is easy to see that $N^{F}(x;\Omega)\subset N^{M}(x;\Omega)$.

######  Proposition 1.

Let $\Omega\subset\mathcal{H}$ be a nonempty closed set, and let
$\delta_{\Omega}$ be the indicator function on $\Omega$. Then we have

$\displaystyle\partial^{M}\delta_{\Omega}(x)=N^{M}(x;\Omega).$

###### Proof.

See [15, Proposition 1.79]. ∎

The following theorem is the fundamental calculus rule.

######  Theorem 1.

Let $f,g:\mathcal{H}\to\mathbb{R}$ be functions such that $f$ is locally
Lipschitz around some point $x$ and $g$ is lower semicontinuous around this
point. Then

$\displaystyle\partial^{M}(f+g)(x)\subset\partial^{M}f(x)+\partial^{M}g(x).$

###### Proof.

See [15, Theorem 3.36]. ∎

The next result ensures that the graph of Mordukhovich’s subdifferential is a
closed set in the norm-weak topology of ${\mathcal{H}}\times{\mathcal{H}}$
when $f$ is locally Lipschitz around the point in question. This fact is
essential in the convergence analysis of the designed algorithm to find Pareto
critical points.

######  Theorem 2.

Let $\\{x^{k}\\},\\{\omega^{k}\\}\subset\mathcal{H}$ be sequences such that
$\omega^{k}\in\partial^{M}f(x^{k})$ for all
$k\in\mathbb{N},\,x^{k}\rightharpoonup x$, and
$\omega^{k}\rightharpoonup\omega$ as $k\to\infty$. If
$f:\mathcal{H}\to\mathbb{R}$ is locally Lipschitz around $x$, then we have the
inclusion $\omega\in\partial^{M}f(x)$.

###### Proof.

It is well known in geometric theory of Banach spaces that any Hilbert space
is Asplund and weakly compactly generated (WCG); see, e.g., [15, Chapter 2]
with the further references and commentaries therein. Thus we are in the
setting of [15, Theorem 3.60]. This theorem tells us that we have the claimed
subdifferential graph closedness if the graph of the subdifferential mapping
$\partial^{M}f$ in $\mathcal{H}\times\mathcal{H}$ enjoys the sequential normal
compactness (SNC) property from [15, Definition 1.20]. This property is fairly
general and holds, in particular, when the set in question is epi-Lipschitzian
in the sense of Rockafellar [16]. As shown in [16], the latter is surely the
case for the graphs of locally Lipschitz functions. Thus we are done with the
proof of the closed-graph theorem. ∎

######  Remark 2.

The aforementioned result of [15, Theorem 3.60] allows us to get the closed-
graph property of $\partial^{M}f$ for significantly larger classes of
functions than just the locally Lipschitz ones. Observe that the SNC property
holds automatically in finite-dimensional spaces, which allows us to cover in
the corresponding counterpart of Theorem 2 a very broad class of extended-
real-valued subdifferentially continuous [17] functions in finite dimensions.
In infinite-dimensional spaces, even appealing to the epi-Lipschitzian
property of the limiting subdifferential graph, we cover the class of
(locally) directionally Lipschitzian functions, which are much broader than
the standard locally Lipschitz class; see [16]. Moreover, we know from [15,
Theorem 1.26] that, in any Banach space, a set is SNC around a given point if
it is compactly epi-Lipschitzian (CEL) around this point in the sense of
Borwein-Strójwas [5]. This creates the opportunity of extending the main
algorithmic and related results of this paper to truly non-Lipschitzian
functions in finite and infinite dimensions.

To proceed with the study of weak Pareto points, we present the following
lemma useful in deriving necessary conditions for weak Pareto optimality.

######  Lemma 1.

Let $\mathcal{D}:=\\{x\ \in\mathcal{C}:g_{j}(x)\leq 0,j\in\mathcal{J}\\}$,
where $g_{j}:\mathcal{H}\to\mathbb{R}$ and $\mathcal{J}=\\{1,\cdots,p\\}$. If
$x^{*}$ is a weak Pareto point of $F$ in $\mathcal{D}$, then for all
$\varepsilon>0$ and $x\in\mathcal{C}$ we have the positivity condition

$\displaystyle
F_{\varepsilon}(x):=\max\\{f_{i}(x)-f_{i}(x^{*})+\varepsilon,g_{j}(x):i\in\mathcal{I},\,j\in\mathcal{J}\\}>0.$

###### Proof.

Suppose the contrary; i.e., there exist $\varepsilon_{0}$ and
$x_{1}\in\mathcal{C}$ such that

$\displaystyle F_{\varepsilon_{0}}(x_{1})\leq 0.$

Therefore, it follows that

$\displaystyle f_{i}(x_{1})<f_{i}(x_{1})+\varepsilon_{0}\leq
f_{i}(x^{*}),\quad i\in\mathcal{I},$ $\displaystyle g_{j}(x_{1})\leq 0,\quad
j\in\mathcal{J},$

and so $x_{1}\in\mathcal{D}$. This contradicts the weak Pareto optimality of
$x^{*}$ and hence verifies the claim of the lemma. ∎

The next result establishes new optimality conditions of the Fritz-John type,
provided by Mordukhovich’s subdifferential, which are necessary for weak
Pareto optimal solutions used in the PPA while being of their own interest.

######  Theorem 3.

Let $\mathcal{C}\subset\mathcal{H}$ be a nonempty, closed, and convex set.
Assume that the set $\mathcal{D}$ is represented in the form

$\mathcal{D}=\\{x\in\mathcal{C}:g_{j}(x)\leq 0,j\in J\\},$

where the functions $f_{i},g_{j}:\mathcal{H}\to\mathbb{R}$ are locally
Lipschitz around the reference point for all $i\in\mathcal{I}$ and
$j\in\mathcal{J}$. If $x^{*}\in\mathcal{D}$ is a weak Pareto solution of
(2.1), then there exist real numbers $\tau>0$ and $\alpha_{i}\geq 0$,
$\beta_{j}\geq 0$ with $i\in\mathcal{I}$ and $j\in\mathcal{J}$ such that we
have the conditions

$0\in\sum_{i\in\mathcal{I}}\alpha_{i}\partial^{M}f_{i}(x^{*})+\sum_{j\in\mathcal{J}}\beta_{j}\partial^{M}g_{j}(x^{*})+\tau\partial^{M}d_{\mathcal{C}}(x^{*}),$
$\mbox{with}\quad\sum_{i\in\mathcal{I}}\alpha_{i}+\sum_{j\in\mathcal{J}}\beta_{j}=1\quad\mbox{and}\quad\beta_{j}g_{j}(x^{*})=0,\quad
j\in\mathcal{J}.$

###### Proof.

Let $x^{*}$ be a weak Pareto point of $F$ in $\mathcal{D}$, and let
$\\{\varepsilon_{l}\\}\subset\mathbb{R}_{++}$ be a sequence converging to 0 as
$l\to\infty$. It follows from Lemma 1 that

$\displaystyle F_{\varepsilon_{l}}(x)>0\;\mbox{ for all
}\;l\in\mathbb{N},\;\;x\in\mathcal{C}.$

Therefore, the function $F_{\varepsilon_{l}}:\mathcal{C}\to\mathbb{R}_{+}$ is
bounded from below; in particular, when $0<\varepsilon_{l}\leq 1$ in our
consideration. The imposed local Lipschitz continuity of $f_{i},g_{j}$ around
$x^{*}$ for all $i\in\mathcal{I}$ and $j\in\mathcal{J}$ ensures that
$F_{\varepsilon_{l}}$ is also locally Lipschitz around this point.
Furthermore, as $x^{*}\in\mathcal{D}$ and $F_{\varepsilon_{l}}(x)>0$ for all
$x\in\mathcal{C}$, we get that
$\displaystyle\inf_{x\in\mathcal{C}}F_{\varepsilon_{l}}(x)\geq 0$, which tells
us in turn that

$\displaystyle\varepsilon_{l}=F_{\varepsilon_{l}}(x^{*})\leq\inf_{x\in\mathcal{C}}F_{\varepsilon_{l}}(x)+\varepsilon_{l}.$

The latter means that $F_{\varepsilon_{l}}$ satisfies all the hypotheses of
[14, Lemma 2.2 ], which is an appropriate version of the fundamental Ekeland
variational principle; see, e.g., [15, Theorem 2.26]. In this way, we find
$z^{l}\in\mathcal{C}$ such that

(2.3) $\displaystyle\|x^{*}-z^{l}\|\leq\sqrt{\varepsilon_{l}}\;\mbox{ and}$
$\displaystyle F_{\varepsilon_{l}}(x)+\sqrt{\varepsilon_{l}}\|x-z^{l}\|\geq
F_{\varepsilon_{l}}(z^{l})\;\mbox{ for all }\;x\in\mathcal{C},$

i.e., $z^{l}$ is a global solution of
$F_{\varepsilon_{l}}(\cdot)+\sqrt{\varepsilon_{l}}\|\cdot-z^{l}\|$ in
$\mathcal{C}$. Denoting by $L$ a Lipschitz constant of
$F_{\varepsilon_{l}}(\cdot)$, we get that $L+1$ is a Lipschitz constant of
$F_{\varepsilon_{l}}(\cdot)+\sqrt{\varepsilon_{l}}\|\cdot-z^{l}\|$. Thus, by
[14, Lemma 2.1], $z^{l}$ is a local optimal solution of the unconstrained
problem of minimizing the function

$\displaystyle
G_{\varepsilon_{l}}(x)=F_{\varepsilon_{l}}(x)+\sqrt{\varepsilon_{l}}\|x-z^{l}\|+(L+1)d_{\mathcal{C}}(x).$

Applying to the latter problem the generalized Fermat rule in terms of the
limiting subdifferential (see [15, Proposition 1.114]) and then the
subdifferential sum rule from Lemma 1, we arrive at

$\displaystyle
0\in\partial^{M}G_{\varepsilon_{l}}(z^{l})\subset\partial^{M}F_{\varepsilon_{l}}(z^{l})+\sqrt{\varepsilon_{l}}\partial^{M}(\|\cdot-z^{l}\|)(z^{l})+(L+1)\partial^{M}d_{\mathcal{C}}(z^{l}).$

Therefore, for each $l$ there exists
$\omega^{l}\in\partial^{M}F_{\varepsilon_{l}}(z^{l})$ satisfying

(2.4) $\displaystyle
0\in\omega^{l}+\sqrt{\varepsilon_{l}}\partial^{M}(\|\cdot-z^{l}\|)(z^{l})+(L+1)\partial^{M}d_{\mathcal{C}}(z^{l}).$

Using [15, Theorem 3.46] and then [15, Theorem 3.52] gives us nonnegative
numbers $u_{i}^{l}$ and $v_{j}^{l}$ such that

(2.5)
$\displaystyle\sum_{i\in\mathcal{I}}\alpha^{l}_{i}+\sum_{j\in\mathcal{J}}\beta^{l}_{j}=1,$
(2.6) $\displaystyle\omega^{l}\in\sum_{\\{i|\alpha^{l}_{i}\neq 0\\}}$
$\displaystyle\alpha_{i}^{l}\partial^{M}f_{i}(z^{l})+\sum_{\\{j|\beta^{l}_{j}\neq
0\\}}\beta_{j}^{l}\partial^{M}g_{j}(z^{l}).$

Combining (2.4) and (2.6) ensures that

$\begin{array}[]{ll}0\in\displaystyle\sum_{\\{i|\alpha^{l}_{i}\neq
0\\}}\alpha_{i}^{l}\partial^{M}f_{i}(z^{l})+\sum_{\\{j|\beta^{l}_{j}\neq
0\\}}\beta_{j}^{l}\partial^{M}g_{j}(z^{l})\\\
\qquad\displaystyle+\sqrt{\varepsilon_{l}}\partial^{M}(\|\cdot-z^{l}\|)+(L+1)\partial^{M}d_{\mathcal{C}}(z^{l}).\end{array}$

In this way, we find sequences
$\\{u_{i}^{l}\\},\,\\{v_{j}^{l}\\},\,\\{\tilde{u}^{l}\\}$, and
$\\{\tilde{v}^{l}\\}$ with

$u_{i}^{l}\in\partial^{M}f_{i}(z^{l}),\;v_{j}^{l}\in\partial^{M}g_{j}(z^{l}),\;\tilde{u}^{l}\in\sqrt{\varepsilon_{l}}\partial^{M}(\|\cdot-z^{l}\|),$

and $\tilde{v}^{l}\in(L+1)\partial^{M}d_{\mathcal{C}}(z^{l})$ such that

(2.7) $\displaystyle 0=\sum_{\\{i|\alpha^{l}_{i}\neq 0\\}}$
$\displaystyle\alpha_{i}^{l}u_{i}^{l}+\sum_{\\{j|\beta^{l}_{j}\neq
0\\}}\beta_{j}^{l}v_{j}^{l}+\tilde{u}^{l}+\tilde{v}^{l}.$

If follows from (2.5) and (2.3) that the sequences $\\{\alpha_{i}^{l}\\}_{l}$,
$\\{\beta_{j}^{l}\\}$, and $\\{z^{l}\\}$ are bounded. Furthermore, the
subgradient boundedness of $\partial^{M}$ by the local Lipschitz continuity of
$f_{i}$, $g_{j}$, and $\|\cdot-z^{l}\|$ (see [15, Corollary 1.81]) ensures the
boundedness of the sequences $\\{u_{i}^{l}\\}$, $\\{v_{j}^{l}\\}$, and
$\\{\tilde{u}^{l}\\}$. The boundedness of $\\{\tilde{v}^{l}\\}$ follows from
(2.7). Passing to subsequences if necessary, we get that the sequences
$\\{\alpha_{i}^{l}\\}$, $\\{\beta_{j}^{l}\\}$, $\\{u_{i}^{l}\\}$,
$\\{v_{j}^{l}\\}$, $\\{\tilde{u}^{l}\\}$, and $\\{\tilde{v}^{l}\\}$ weakly
converge, respectively, to $\alpha_{i}$, $\beta_{j}$, $u_{i}$, $v_{j}$,
$\tilde{u}$, and $\tilde{v}$. Hence $\alpha_{i},\;\beta_{j}\geq 0$ and
$\sum_{i\in\mathcal{I}}\alpha_{i}+\sum_{j\in\mathcal{J}}\beta_{j}=1$. By
(2.3), we have that $z^{l}\to x^{*}$. Then (2.7) as $l\to\infty$ gives us

$\displaystyle
0=\sum_{i\in\mathcal{I}}\alpha_{i}u_{i}+\sum_{j\in\mathcal{J}}\beta_{j}v_{j}+\tilde{u}+\tilde{v}.$

The closed-graph property of Theorem 2 yields

$\displaystyle u_{i}\in\partial^{M}f_{i}(x^{*}),\quad
v_{j}\in\partial^{M}g_{j}(x^{*}),\quad\mbox{and}\quad\tilde{v}\in\partial^{M}d_{\mathcal{C}}(x^{*}),$

and therefore, we arrive at the claimed inclusion

$\displaystyle
0\in\sum_{i\in\mathcal{I}}\alpha_{i}\partial^{M}f_{i}(x^{*})+\sum_{j\in\mathcal{J}}\beta_{j}\partial^{M}g_{j}(x^{*})+\tau\partial^{M}d_{\mathcal{C}}(x^{*})\;\mbox{
with some }\;\tau>0.$

∎

To proceed further, recall the following result from [15, Theorem 1.97].

######  Lemma 2.

Let $\Omega\subset\mathcal{H}$ be nonempty and closed. Then

$N^{M}(x;\Omega)=\displaystyle\cup_{\lambda>0}\lambda\partial^{M}d_{\Omega}(x)\;\mbox{
for any }\;x\in\Omega.$

Here is the subdifferential estimate for the distance function, which plays an
important role in the subsequent convergence analysis.

######  Theorem 4.

Let $\Omega$ be a nonempty, closed, and convex subset of $\mathcal{H}$. Then
for any $x\in\Omega$, we have the upper estimate

$\displaystyle\partial^{M}d_{\Omega}(x)\subset B[0,1]\cap N^{M}(x;\Omega),$

where $B[0,1]$ denotes the closed unit ball in $\mathcal{H}$.

###### Proof.

It is known that the distance function $d_{\Omega}(x)$ is globally Lipschitz
with a constant $L=1$. It follows from [15, Corollary 1.81] that
$\partial^{M}d_{\Omega}(x)\subset B[0,1]$. Employing now Lemma 2 tells us that
$N^{M}(x;\Omega)=\displaystyle\cup_{\lambda>0}\lambda\partial^{M}d_{\Omega}(x)$
, and thus we arrive at $\partial^{M}d_{\Omega}(x)\subset N^{M}(x;\Omega)$. ∎

##  3\.  The proximal point algorithm

In this section, we prove some facts related to our new approach for
convergence of the Proximal Point Algorithm (PPA) for vector-valued functions.

###  3.1. The algorithm

In what follows, $\Omega\subset\mathcal{H}$ is a nonempty, closed, and convex
set while $F=(f_{1},\ldots,f_{m}):\mathcal{H}\to\mathbb{R}^{m}$ is a vector
function such that its each component $f_{i}:H\to\mathbb{R}$,
$i\in\mathcal{I}$, is locally Lipschitz. From Remark 1, we can assume without
loss of generality that $F\succ 0$.

To formulate the PPA for finding a Pareto critical point of $F$ in $\Omega$,
let $\\{\lambda_{k}\\}$ be a sequence of positive numbers, and let
$\\{\varepsilon^{k}\\}\subset\mathbb{R}_{++}^{m}$ be a sequence with
$||\varepsilon^{k}||=1$ for all $k\in\mathbb{N}$. The method generates an
iterative sequence $\\{x^{k}\\}\subset\Omega$ in the following way:

Algorithm (refined vectorial PPA)
Initialization: Choose $x^{0}\in\Omega$.
Stopping rule: If for $x^{k}\in\Omega$ we have that $x^{k}$ is a Pareto
critical point, then set $x^{k+p}:=x^{k}$ for all $p\in\mathbb{N}$.
Iterative step: As the next iterate, take $x^{k+1}\in\Omega$ such that

$x^{k+1}\in\mbox{argmin}_{w}\left\\{F(x)+\frac{\lambda_{k}}{2}||x-x^{k}||^{2}\varepsilon^{k}\;:\;x\in\Omega_{k}\right\\},$

where $\Omega_{k}:=\\{x\in\Omega\;:\;F(x)\preceq F(x^{k})\\}$.

We would like to mention that this method aims at finding separate solutions
while not the entire solution set. It has been noticed by Fukuda and Graña
Drummond [10] and by Fliege et al. [9] that we can expect to somehow
approximate the solution set by just performing this method for different
initial points. In the well-known weighting method, this kind of idea also
appears in Burachik et al. [6]. More precisely, the method can be performed
for different weights in order to find either the solution set, or a
reasonable approximation of this set. However, in some cases, arbitrary
choices of the weighting vectors may lead the weighting method to unbounded
problems. The Pareto front, i.e., the objective values of these solutions, is
in general an infinite set. Hence, in practice, only an approximation of the
Pareto front is obtained.

######  Proposition 2.

The Algorithm is well-defined.

###### Proof.

The starting point $x^{0}\in\Omega$ is chosen in the initialization step.
Assuming in iteration $k$ that $x^{k}$ is given and denote

$F_{k}(x):=F(x)+\frac{\lambda_{k}}{2}\|x-x^{k}\|^{2}\varepsilon^{k}.$

Note that we have $x^{k}\in\Omega_{k}$, which implies that $F_{k}(\Omega_{k})$
is nonempty. It is straightforward to check that $F_{k}(x)\succ 0$ and
$F_{k}(\Omega_{k})$ is closed. Since the cone $\mathbb{R}^{m}_{+}$ is Daniell,
it follows from [13, Lemma 3.5] that $F_{k}(\Omega_{k})$ is
$\mathbb{R}^{m}_{+}$-complete. From [13, Theorem 3.3], the next iteration
exists and satisfies the inclusion

$x^{k+1}\in\mbox{arg\,min}_{w}\\{F_{k}(x)\,:\,x\in\Omega_{k}\\}.$

which completes the proof of the proposition. ∎

From now on, the sequences $\\{x^{k}\\}$, $\\{\lambda_{k}\\}$, and
$\\{\varepsilon^{k}\\}$ are taken from the above Algorithm. The following
hold.

######  Proposition 3.

For all $k\in\mathbb{N}$, there exist
$\alpha^{k},\beta^{k}\in\mathbb{R}^{m}_{+}$,
$\omega^{k},u_{i}^{k}\in\mathcal{H}$, and $\tau_{k}\in\mathbb{R}_{++}$
satisfying the equalities

(3.1)
$\sum_{i\in\mathcal{I}}(\alpha^{k}_{i}+\beta^{k}_{i})u_{i}^{k}+\lambda_{k-1}\sum_{i\in\mathcal{I}}\varepsilon^{k-1}_{i}\alpha^{k}_{i}(x^{k}-x^{k-1})+\tau_{k}\omega^{k}=0$

together with the conditions

$\quad\omega^{k}\in B[0,1]\cap N^{M}(x^{k};\Omega),\quad
u^{k}_{i}\in\partial^{M}f_{i}(x^{k}),\quad\mbox{and}\quad||\alpha^{k}+\beta^{k}||_{1}=1.$

###### Proof.

It follows from the constructions of the Algorithm that $x^{k}$ is a weak
Pareto solution of the multiobjective problem

$\min_{w}\\{F_{k-1}(x)\;:\;x\in\Omega_{k-1}\\}$,

where
$F_{k-1}(x):=F(x)+\frac{\lambda_{k-1}}{2}||x-x^{k-1}||^{2}\varepsilon^{k-1}$
for each $k\in\mathbb{N}$.

Defining $G_{k-1}(x):=F(x)-F(x^{k-1})$, it is easy to deduce from the local
Lipschitz continuity of $F$ that all the component functions

(3.2) $(g_{k-1})_{i}(\cdot)=f_{i}(\cdot)-f_{i}(x^{k-1}),\qquad
i\in\mathcal{I},$

with $(f_{k-1})_{i}$ given by

(3.3)
$(f_{k-1})_{i}(\cdot)=f_{i}(\cdot)+\frac{\lambda_{k-1}}{2}||\cdot-x^{k-1}||^{2}\varepsilon^{k-1}_{i},\qquad
i\in\mathcal{I},$

are locally Lipschitzian. Fix $k\in\mathbb{N}$ and apply Lemma 1 and Theorem 3
with $g_{j}$ and $f_{j}$ defined in $\eqref{gtilde}$ and $\eqref{ftilde}$,
respectively. Taking into account that $\partial^{M}d_{\Omega}(x^{k})\subset
B[0,1]\cap N^{M}(x^{k};\Omega)$, $k\in\mathbb{N}$, by Theorem 4 ensures that

$\displaystyle
0=\sum_{i\in\mathcal{I}}(\alpha_{i}^{k}+\beta_{i}^{k})u_{i}^{k}+\lambda_{k-1}\sum_{i\in\mathcal{I}}\alpha_{i}^{k}\varepsilon_{i}^{k-1}(x^{k}-x^{k-1})+\tau_{k}\omega^{k},$

where the elements involved satisfy the conditions

$u_{i}^{k}\in\partial^{M}f_{i}(x^{k}),\quad\omega^{k}\in\partial^{M}d_{\Omega}(x^{k})\subset
B[0,1]\cap
N^{M}(x^{k};\Omega),\quad\displaystyle\sum_{i\in\mathcal{I}}(\alpha_{i}^{k}+\beta_{i}^{k})=1$

for all $k$. This completes the proof of the proposition. ∎

As in [4], the stopping rule in the Algorithm can be changed by the following
rule, which is easier to check: after computing $x^{k+1}$, the Algorithm stops
if $x^{k+1}=x^{k}$, i.e., we set $x^{k+p}:=x^{k}$ for all $p\geq 1$.

######  Corollary 1.

If $x^{k+1}=x^{k}$, then $x^{k}$ is a Pareto critical point of $F$.

###### Proof.

Assume that $x^{k+1}=x^{k}$. Then Proposition 3 tells us that there exist
$\alpha^{k+1},\;\beta^{k+1}\in\mathbb{R}^{m}_{+},\;u_{i}^{k+1},\omega^{k+1}\in\mathcal{H}$,
and $\tau_{k+1}\in\mathbb{R}_{++}$ satisfying

$\displaystyle\sum_{i\in\mathcal{I}}(\alpha_{i}^{k+1}+\beta_{i}^{k+1})u_{i}^{k+1}+\tau_{k+1}\omega^{k+1}=0,$

where $-\sum_{i\in\mathcal{I}}(\alpha_{i}^{k+1}+\beta_{i}^{k+1})u_{i}^{k+1}\in
N^{M}(x^{k+1};\Omega)$. If $x^{k}$ is not a Pareto critical point of $F$ in
$\Omega$, then there exists $y\in\Omega$ such that $\langle
u_{i}^{k+1},y-x^{k+1}\rangle<0$. Therefore, with
$\displaystyle\sum_{i\in\mathcal{I}}(\alpha_{i}^{k+1}+\beta_{i}^{k+1})=1$ we
get

$\displaystyle\left\langle\sum_{i\in\mathcal{I}}(\alpha_{i}^{k+1}+\beta_{i}^{k+1})u_{i}^{k+1},y-x^{k+1}\right\rangle<0,$

which contradicts
$-\sum_{i\in\mathcal{I}}(\alpha_{i}^{k+1}+\beta_{i}^{k+1})u_{i}^{k+1}\in
N^{M}(x^{k+1};\Omega)$ and thus verifies that $x^{k}$ is a Pareto critical
point of $F$ in $\Omega$. ∎

######  Corollary 2.

If there exists a number $k_{0}\in\mathbb{N}$ such that $\alpha^{k_{0}}=0$,
then $x^{k_{0}}$ is a Pareto critical point of $F$.

###### Proof.

Assume that there exists $k_{0}\in\mathbb{N}$ such that $\alpha^{k_{0}}=0$.
Then Proposition 3 provides the equality

$\displaystyle\sum_{i\in\mathcal{I}}\beta_{i}^{k_{0}}u_{i}^{k_{0}}+\tau^{k_{0}}\omega^{k_{0}}=0$

with $\omega^{k_{0}}\in B[0,1]\cap N^{M}(x^{k_{0}};\Omega)$,
$u_{i}^{k_{0}}\in\partial^{M}f_{i}(x^{k_{0}})$, and
$\displaystyle\sum_{i\in\mathcal{I}}\beta_{i}^{k_{0}}=1$. Therefore, we get
$\beta^{k_{0}}\in\mathbb{R}^{m}_{+}$. Supposing now that $x^{k_{0}}$ is not a
Pareto critical point of $F$ in $\Omega$ ensures the existence of $y\in\Omega$
such that for all $i\in\mathcal{I}$ and
$u_{i}^{k_{0}}\in\partial^{M}f_{i}(x^{k_{0}})$ we have the strict inequality

$\displaystyle\langle u_{i}^{k_{0}},y-x^{k_{0}}\rangle<0.$

Finally, it follows from $\beta_{i}^{k_{0}}\in\mathbb{R}^{m}_{+}$ and
$\displaystyle\sum_{i\in\mathcal{I}}\beta_{i}^{k_{0}}=1$ that

$\displaystyle\left\langle\sum_{i\in\mathcal{I}}\beta_{i}^{k_{0}}u_{i}^{k_{0}},y-x^{k_{0}}\right\rangle<0,$

which contradicts the fact that
$-\displaystyle\sum_{i\in\mathcal{I}}\beta_{i}^{k_{0}}u_{i}^{k_{0}}\in
N^{M}(x^{k_{0}};\Omega)$. ∎

######  Remark 3.

Note that if the Algorithm terminates after a finite number of iterations, it
terminates at a Pareto critical point. This allows us us to suppose that
$\\{x^{k}\\}$ is an infinite sequence, and hence $x^{k+1}\neq x^{k}$ for all
$k\in\mathbb{N}$ due to the result of Corollary 1.

##  4\. Convergence analysis

In this section, we employ the following $\mathbb{R}^{m}_{+}$-completeness
assumption imposed on the vector function $F$ and the initial point $x^{0}$:

($\mathbb{R}^{m}_{+}$-Completeness assumption) For all sequences
$\\{a^{k}\\}\subset\mathcal{H}$ satisfying $a^{0}=x^{0}$ and
$F(a^{k+1})\preceq F(a^{k})$ for all $k\in\mathbb{N}$, there exists
$a\in\mathcal{H}$ such that

$F(a)\preceq F(a^{k})\;\mbox{ for all }\;k\in\mathbb{N}.$

In the classical case of unconstrained convex optimization, this condition is
equivalent to the existence of solutions to the optimization problem in
question. In general, the $\mathbb{R}^{m}_{+}$-completeness assumption is
recognized to ensure the existence of efficient points for vector optimization
problems; see, e.g., [13, page 46]. An interesting discussion on the existence
of efficient points for vector optimization problems can also be found in [13,
Chapter 2]. In what follows, we employ this conditions to establishing the
convergence of our Algorithm to find Pareto critical points in problems of
multiobjective optimization.

###  4.1. Locally Lipschitz case

First we establish our convergence result in the locally Lipschitz case.

######  Theorem 5.

Suppose that there exist scalars $a,b,c\in\mathbb{R}_{++}$ such that
$0<a\leq\lambda_{k}\leq b$ and $0<c\leq\varepsilon_{j}^{k}$, for all
$k\in\mathbb{N}$ and $j=1,\ldots,m$. Then every weak cluster point of
$\\{x^{k}\\}$ is a Pareto critical point of $F$.

###### Proof.

It follows from the definition of the Algorithm that $x^{k}$, for each fixed
$k\in\mathbb{N}$, is an optimal solution of the problem

$\displaystyle\mbox{min}_{w}\left\\{F(x)+\frac{\lambda_{k-1}}{2}\|x-x^{k-1}\|^{2}\varepsilon^{k-1}\;:\;x\in\Omega_{k-1}\right\\}.$

This readily yields the condition

$\displaystyle\max_{1\leq j\leq
m}\left\\{f_{j}(x^{k-1})-f_{j}(x^{k})-\frac{\lambda_{k-1}}{2}\|x^{k}-x^{k-1}\|^{2}\varepsilon^{k-1}\right\\}\geq
0.$

Take $j_{0}(k)=j_{0}\in\\{1,...,m\\}$ as an index where the maximum in the
last inequality is attained. Then we deduce from the lower boundedness
assumption on the sequences $\\{\lambda_{k}\\}$ and $\\{\varepsilon^{k}\\}$
that

$\displaystyle\dfrac{ac}{2}\|x^{k}-x^{k-1}\|$ $\displaystyle\leq
f_{j_{0}}(x^{k-1})-f_{j_{0}}(x^{k})$
$\displaystyle=\sqrt{(f_{j_{0}}(x^{k-1})-f_{j_{0}}(x^{k}))^{2}}$
$\displaystyle\leq\sqrt{\sum_{j\in\mathcal{I}}(f_{j}(x^{k-1})-f_{j}(x^{k}))^{2}}$
(4.1) $\displaystyle=\|F(x^{k-1})-F(x^{k})\|,\quad k\in\mathbb{N}.$

Since $\\{F(x^{k})\\}$ is nonincreasing and $F\succ 0$, we deduce that the
right-hand side of (4.1) converges to $0$ as $k\to\infty$. This gives us

(4.2) $\displaystyle\|x^{k}-x^{k-1}\|\to 0\;\;as\;\;k\to\infty.$

Now let $\hat{x}$ be a weak cluster point of $\\{x^{k}\\}$, and let
$\\{x^{k_{l}}\\}$ be a subsequence of $\\{x^{k}\\}$ converging to $\hat{x}$.
Applying Proposition 2 to the sequence $\\{x^{k_{l}}\\}$, we find

$\left\\{\begin{array}[]{l}\\{u_{i}^{k_{l}}\\}\subset\partial^{M}f_{i}(x^{k_{l}}),\\\
\\{\alpha^{k_{l}}\\},\\{\beta^{k_{l}}\\}\subset\mathbb{R}^{m}_{+},\\\
\\{\omega^{k_{l}}\\}\subset\mathcal{H},\quad\mbox{and}\quad\\{\tau_{k_{l}}\\}\subset\mathbb{R}_{++}\end{array}\right.$

satisfying the equalities

(4.3)
$\displaystyle\sum_{i\in\mathcal{I}}(\alpha_{i}^{k_{l}}u_{i}^{k_{l}}+\beta_{i}^{k_{l}})+\lambda_{k_{l}-1}\sum_{i\in\mathcal{I}}\varepsilon^{k_{l}-1}_{i}\alpha^{k_{l}}_{i}(x^{k_{l}}-x^{k_{l}-1})+\tau_{k_{l}}\omega^{k_{l}}=0,$

for all $k\in\mathbb{N}$, where

(4.4) $\displaystyle\omega^{k_{l}}\in B[0,1]\cap
N^{M}(x^{k_{l}};\Omega)\quad\mbox{and}\quad\displaystyle\sum_{i\in\mathcal{I}}(\alpha^{k_{l}}+\beta^{k_{l}})=1.$

Since $\\{x^{k_{l}}\\}$ weakly converges to $\hat{x}$, it follows that
$\\{x^{k_{l}}\\}$ is bounded. Using the local Lipschitz continuity of $F$ and
taking into account that $\partial^{M}f_{j}$ in (4.4) is bounded on compact
sets by [15, Corollary 1.81], we get that the sequences
$\\{\alpha^{k_{l}}\\},\,\\{\beta^{k_{l}}\\},\,\\{\omega^{k_{l}}\\}$ are
bounded. Thus $\\{\tau_{k_{l}}\\}$ is bounded as well. Without loss of
generality, assume that
$\\{u_{i}^{k_{l}}\\},\,\\{\alpha^{k_{l}}\\},\,\\{\beta^{k_{l}}\\},\,\\{\omega^{k_{l}}\\}$
weakly converge to $u_{i},\,\alpha,\,\beta,\,\omega$, respectively, and that
$\\{\tau_{k_{l}}\\}$ converge to $\tau$. Furthermore, it follows from the
boundedness of
$\lambda_{k_{l}-1}\sum_{i\in\mathcal{I}}\varepsilon^{k_{l}-1}_{i}\alpha^{k_{l}}_{i}$
together with (4.2) that

$\displaystyle\sum_{i\in\mathcal{I}}(\alpha_{i}+\beta_{i})u_{i}+\tau\omega=0.$

Employing Lemma 1 and Theorem 2 gives us $u_{i}\in\partial^{M}f_{i}(\hat{x})$
and $\omega\in N^{M}(\hat{x};\Omega)$. Thus we arrive at the inclusion

(4.5) $\displaystyle-\sum_{i\in\mathcal{I}}(\alpha_{i}+\beta_{i})u_{i}\in
N^{M}(\hat{x};\Omega).$

Arguing finally by contradiction, suppose that $\hat{x}$ is not a Pareto
critical point of $F$. Then there exists $y\in\Omega$ such that

$\displaystyle\langle u_{i},y-\hat{x}\rangle<0\;\mbox{ for all
}\;i\in\mathcal{I}.$

This implies therefore the strict inequality

$\left\langle\displaystyle\sum_{i\in\mathcal{I}}(\alpha_{i}+\beta_{i})u_{i},y-\hat{x}\right\rangle<0.$

Remembering that $\alpha,\beta\in\mathbb{R}^{m}_{+}$ and
$\displaystyle\sum_{i\in\mathcal{I}}(\alpha_{i}+\beta_{i})=1$, we get a
contradiction with (4.5) and hence complete the proof of the theorem. ∎

###  4.2. Pseudoconvexity case

This subsection presents a constructive sufficient condition ensuring that the
entire sequence of the iterates generated by the Algorithm weakly converges to
a Pareto critical point of $F$.

Using the limiting subdifferential $\partial^{M}$, we say that a scalar
function $f$ defined on a Hilbert space ${\mathcal{H}}$ is pseudoconvex if for
any $x,y\in\mathcal{H}$ we have

$\displaystyle\exists\,\omega\in\partial^{M}f(x):\langle\omega,y-x\rangle\geq
0\Rightarrow f(x)\leq f(y).$

Similarly to the differentiable case, it is not difficult to observe that any
pseudoconvex function $f$ satisfies the following fundamental properties:

  * (a)

every local minimizer of $f$ is a global minimizer;

  * (b)

$0\in\partial^{M}f(x)\Rightarrow f$ has a global minimum at $x$.

Further, we say that the vector function
$F=(f_{1},\ldots,f_{m}):\mathcal{H}\to\mathbb{R}^{m}$ is pseudoconvex if each
component function $f_{i}:H\to\mathbb{R}$, $i\in\mathcal{I}$, is a
pseudoconvex.

######  Proposition 4.

Let $F:\mathcal{H}\to\mathbb{R}^{m}$ be a pseudoconvex function. If
$x\in\Omega$ is a Pareto critical point of $F$ in $\Omega$, then it is a weak
Pareto solution of $F$ in $\Omega$.

###### Proof.

Assume that $x$ is a Pareto critical point of $F$ in $\Omega$. Arguing by
contradiction, suppose that $x$ is not a weak Pareto solution of $F$ in
$\Omega$. Then there exists a point $\tilde{x}\in\Omega$ such that

$\displaystyle f_{i}(\tilde{x})<f_{i}(x)\;\mbox{ whenever }\;i\in\mathcal{I}.$

Since $f_{i}$ is pseudoconvex, it follows that for all
$\omega_{i}\in\partial^{M}f_{i}(x)$ and $i\in\mathcal{I}$ we have the strict
inequality

$\displaystyle\langle\omega_{i},\tilde{x}-x\rangle<0.$

However, this contradicts the fact that $x$ is a Pareto critical point of $F$
in $\Omega$. ∎

Before deriving the main result of this subsection, recall that a sequence
$\\{y^{k}\\}$ is said to be Fejér convergent (or Fejér monotone) to a nonempty
set $U\subset\mathcal{H}$ if for all $k\in\mathbb{N}$ we have

$||y^{k+1}-y||\leq||y^{k}-y||\;\mbox{ whenever }\;y\in U.$

The following result is well known, and its proof is elementary.

######  Proposition 5.

Let $U\subset\mathcal{H}$ be a nonempty set, and let $\\{y^{k}\\}$ be a Fejér
convergent sequence to $U$. Then $\\{y^{k}\\}$ is bounded. Moreover, if a weak
cluster point $y$ of $\\{y^{k}\\}$ belongs to $U$, then $\\{y^{k}\\}$
converges to $y$ as $k\to\infty$.

The next theorem shows that in the pseudoconvex case, we have the weak
convergence of the entire sequence of iterations generated by the Algorithm to
a Pareto critical point.

######  Theorem 6.

Under the imposed assumptions, the sequence $\\{x^{k}\\}$ generated by the
Algorithm weakly converges to a weak Pareto solution of F in $\Omega$ as
$k\to\infty$.

###### Proof.

We split the proof of the theorem into the following five steps.
Step 1 (Fejér convergence) Define the set $E\subset\Omega$ by

$E:=\\{x\in{\Omega}\;:\;F(x)\preceq F(x^{k})\;\mbox{ for all
}\;k\in\mathbb{N}\\}.$

From the $\mathbb{R}^{m}_{+}$-completeness assumption on $F$ at $x^{0}$, the
set $E$ is nonempty. Take an arbitrary point $x^{*}\in E$, which means that
$x^{*}\in\Omega_{k}$, and for all $k\in\mathbb{N}$ denote
$\gamma_{k+1}:=\lambda_{k}\langle\varepsilon^{k},\alpha^{k+1}\rangle$. Observe
that $\gamma_{k+1}>0$ whenever $k\in\mathbb{N}$ due to $\lambda_{k}>0$,
$\varepsilon^{k}\in\mathbb{R}^{m}_{++}$ and
$\alpha^{k}\in\mathbb{R}^{m}_{+}\backslash\\{0\\}$. Since

$||x^{k}-x^{*}||^{2}=||x^{k}-x^{k+1}||^{2}+||x^{k+1}-x^{*}||^{2}+2\langle
x^{k}-x^{k+1},x^{k+1}-x^{*}\rangle,$

we deduce from (3.1) that

(4.6)
$\begin{array}[]{ll}\|x^{k}-x^{*}\|^{2}=\|x^{k}-x^{k+1}\|^{2}+\|x^{k+1}-x^{*}\|^{2}\\\
+\displaystyle\dfrac{2}{\theta_{k+1}}\sum_{i\in\mathcal{I}}(\alpha_{i}^{k+1}+\displaystyle\beta_{i}^{k+1})\langle
u_{i}^{k+1},x^{k+1}-x^{*}\rangle\\\
+\dfrac{2}{\theta_{k+1}}\tau_{k+1}\langle\omega_{i}^{k+1},x^{k+1}-x^{*}\rangle,\end{array}$

where $u_{i}^{k+1}\in\partial^{M}f_{i}(x^{k+1})$ for all $k\in\mathbb{N}$ and
$i\in\mathcal{I}$ with
$\theta_{k+1}:=\lambda_{k+1}\sum_{i\in\mathcal{I}}\alpha_{i}^{k+1}\varepsilon_{i}^{k}$.
On the other hand, it follows from the pseudoconvexity of $F$ together with
$x^{*}\in\Omega_{k}$ and $\theta_{k}>0$ for all $k$ that we have

(4.7)
$\frac{2}{\theta_{k+1}}\sum_{i=1}^{m}(\alpha_{i}^{k+1}+\beta_{i}^{k+1})\langle
u_{i}^{k+1},x^{k+1}-x^{*}\rangle\geq 0.$

The inclusion $w^{k+1}\in N^{M}(x^{k+1};\Omega)$ combined with $\tau_{k}>0$
leads us to

(4.8) $\tau_{k+1}\langle w^{k+1},x^{k+1}-x^{*}\rangle\geq 0.$

Using now (4.7) and (4.8) in (4.6) yields

$||x^{k+1}-x^{k}||^{2}\leq||x^{k}-x^{*}||^{2}-||x^{k+1}-x^{*}||^{2},\quad
k\in\mathbb{N},$

which implies that $||x^{k+1}-x^{*}||\leq||x^{k}-x^{*}||$ for any $x^{*}\in
E$. This verifies that $\\{x^{k}\\}$ is Fejér convergent to $E$.
Step 2 (weak cluster points of $\\{x^{k}\\}$ belongs to $E$) Since
$\\{x^{k}\\}$ is Fejér convergent to $E$, it follows from Proposition 5 that
$\\{x^{k}\\}$ is bounded. Let $x^{*}$ be a weak cluster point of
$\\{x^{k}\\}$. The construction of the Algorithm tells us that
$F(x^{k+1})\preceq F(x^{k})$ whenever $k\in\mathbb{N}$. Thus it follows from
the continuity of $F$ that $F(x^{*})\preceq F(x^{k})$ for all $k$, which
verifies that $x^{*}\in E$.
Step 3 (Convergence of the iterative sequence) This step directly follows from
Proposition 5 combined with Steps 1 and 2.
Step 4 (Proximity of consecutive iterates) Assume that $\\{x^{k}\\}$ converges
to $\hat{x}$. It follows from the triangle inequality that

(4.9) $||x^{k+1}-x^{k}||\leq||x^{k+1}-\hat{x}||+||x^{k}-\hat{x}||,\quad
k\in\mathbb{N},$

due to $x^{k}\to\hat{x}$ as $k\to\infty$, we conclude that

$\lim_{k\to\infty}||x^{k+1}-x^{k}||=0,$

which justifies the claim of this step.
Step 5 (The limit point is a weak Pareto solution) The proof of this step uses
the same argument as in the proof of Theorem 5 from (4.3) on. This completes
the proof of the theorem. ∎

######  Acknowledgement.

Research of Boris S. Mordukhovich was partly supported by the US National
Science Foundation under grant DMS-2204519, by the Australian Research Council
under Discovery Project DP190100555, and by the Project 111 of China under
grant D21024. Research of Glaydston de C. Bento was supported in part by CNPq
grants 314106/2020-0 and the research of João Xavier da Cruz Neto was partly
supported by CNPq grants 302156/2022-4.

## References

  * [1] H. C. F. APOLINÁRIO, E. A. PAPA QUIROZ, and P. R. OLIVEIRA, A scalarization proximal point method for quasiconvex multiobjective minimization. J. Glob. Optim., 64 (2016), pp. 79-96.
  * [2] G. C. BENTO, J. X. CRUZ NETO and A. SOUBEYRAN, A proximal point-type method for multicriteria optimization. Set-Valued Var. Anal., 22 (2014), pp. 557–573.
  * [3] G. C. BENTO, J. X. CRUZ NETO, G. LÓPEZ, A. SOUBEYRAN, and J. C. O. SOUZA, The proximal point method for locally Lipschitz functions in multiobjective optimization with application to the compromise problem. SIAM J. Optim., 28 (2018), pp. 1104–1120.
  * [4] H. BONNEL, A. N. IUSEM, and B. F. SVAITER, Proximal methods in vector optimization. SIAM J. Optim., 15 (2005), pp. 953-970.
  * [5] J. M. BORWEIN and H. M. STRÓJWAS, Tangential approximations. Nonlinear Anal., 9 (1985), pp. 1347–1366.
  * [6] R. S. BURACHIK, C. Y. KAYA and M. M. RIZVI, A new scalarization technique to approximate Pareto fronts of problems with disconnected feasible sets. J. Optim. Theory Appl., 162 (2014), pp. 428-446.
  * [7] L. C. CENG and J. C. YAO, Approximate proximal methods in vector optimization. Eur. J. Oper. Res., 183 (2007), pp. 1-19.
  * [8] T. D. CHOUNG, T. D., B. S. MORDUKHOVICH and J. C. YAO, Hybrid approximate proximal algorithms for efficient solutions in vector optimization. J. Nonlinear Convex Anal., 12 (2011), pp. 257-286.
  * [9] J. FLIEGE, L. M. GRAÑA DRUMMOND, and B. F. SVAITER, Newton’s method for multiobjective optimization. SIAM J. Optim., 20 (2009), pp. 602-626.
  * [10] E. H. FUKUDA, AND L. M. GRAÑA DRUMMOND, A survey on multiobjective decent methods. Pesqui. Oper., 34 (2014), pp. 585-620.
  * [11] R. GREGÓRIO and P. R. OLIVEIRA, A logarithmic-quadratic proximal point scalarization method for multiobjective programming. J. Glob. Optim., 49 (2011), pp. 281-291.
  * [12] X. X. HUANG and X. Q. YANG, Duality for multiobjective optimization via nonlinear Lagrangian functions. J. Optim. Theory Appl., 120 (2004), pp. 111–127.
  * [13] D. T. LUC, Theory of Vector Optimization, Lecture Notes in Econom. and Math. Syst. Springer, New York (1989).
  * [14] M. MINAMI, Weak Pareto-optimal necessary conditions in a nondifferentiable multiobjective program on a Banach Space. J. Optim. Theory Appl., 41 (1983), pp. 451–461.
  * [15] B. S. MORDUKHOVICH, Variational Analysis and Generalized Differentiation I: Basic Theory. Grundlehren Series [Fundamental Principles of Mathematical Sciences]. Springer, Berlin (2006).
  * [16] R. T. ROCKAFELLAR, Directional Lipschitzian functions and subdifferential calculus. Proc. London Math. Soc., 39 (1979), pp. 331–355.
  * [17] R. T. ROCKAFELLAR and R. J-B. WETS, Variational Analysis. Grundlehren Series [Fundamental Principles of Mathematical Sciences]. Springer, Berlin (1998).
  * [18] K. D. V. VILLACORTA and P. R. OLIVEIRA, An interior proximal method in vector optimization. Eur. J. Oper. Res., 214 (2011), pp. 485–492.


# The optimal division between sample and background measurement time for
photon counting experiments

Brian Richard Pauw<EMAIL_ADDRESS>–  Samuel Tardif<EMAIL_ADDRESS>–

(August 27, 2024)

##  I Outline

_Note by the author: It was found after this derivation, that the work
presented here had been derived differently but with similar conclusions by
Steinhart and Plestil Steinhart and Plestil (1993). We nevertheless believe
that this simplified derivation may be of use to some readers._

Usually, equal time is given to measuring the background and the sample, or
even a longer background measurement is taken as it has so few counts. While
this seems the right thing to do, the relative error after background
subtraction improves when more time is spent counting the measurement with the
highest amount of scattering. As the available measurement time is always
limited, a good division must be found between measuring the background and
sample, so that the uncertainty of the background-subtracted intensity is as
low as possible.

Herein outlined is the method to determine how best to divide measurement time
between a sample and the background, in order to minimize the relative
uncertainty. Also given is the relative reduction in uncertainty to be gained
from the considered division. It is particularly useful in the case of
scanning diffractometers, including the likes of Bonse-Hart cameras, where the
measurement time division for each point can be optimized depending on the
signal-to-noise ratio.

The optimum setting for machines with photon-counting two-dimensional
detectors has to be further evaluated, but the intention is to include that in
this note.

##  II The calculation

We assume that the number of background intensity photons $I_{b}$, measured
for a time $t_{b}$ is subtracted from the sample measurement photon count
$I_{s}$ which was measured for a time $t_{s}$, to result in the background
subtracted count rate $C_{bs}$:

$C_{bs}=\frac{I_{s}}{t_{s}}-\frac{I_{b}}{t_{b}}$ (1)

Defining the sample uncertainty as $\Delta I_{s}=\sqrt{I_{s}}$ and the
background uncertainty similarly as $\Delta I_{b}=\sqrt{I_{b}}$, the
uncertainty $\Delta C_{bs}$ would then be:

$\Delta C_{bs}=\sqrt{\left(\frac{\Delta
I_{s}}{t_{s}}\right)^{2}+\left(\frac{\Delta
I_{b}}{t_{b}}\right)^{2}}=\sqrt{\frac{I_{s}}{t_{s}^{2}}+\frac{I_{b}}{t_{b}^{2}}}$
(2)

Defining the number of counted photons $I$ to be a multiplication of the count
rate $C$ and the measurement time, we get $I_{b}=C_{b}t_{b}$ and
$I_{s}=C_{s}t_{s}$. Further defining the signal-to-noise ratio
$g=\frac{C_{s}}{C_{b}}$, the total time $t_{t}=t_{b}+t_{s}$ and the fraction
of time spent measuring the sample $f=\frac{t_{s}}{t_{t}}$, we can define our
relative uncertainty in terms of signal-to-noise ratio and time fraction:

$\frac{\Delta
C_{bs}}{C_{bs}}=\sqrt{\frac{\frac{C_{s}}{t_{s}}+\frac{C_{b}}{t_{b}}}{(C_{s}-C_{b})^{2}}}=\sqrt{\frac{1}{C_{b}t_{t}}}\sqrt{\frac{\frac{g}{f}+\frac{1}{(1-f)}}{(g-1)^{2}}}$
(3)

We then can try to find the optimum by locating the value for $f$ where the
derivative of equation LABEL:eq:sErel2 is zero:

$\frac{\partial\frac{\Delta C_{bs}}{C_{bs}}}{\partial
f}=\frac{\partial}{\partial
f}\sqrt{\frac{1}{C_{b}t_{t}}}\sqrt{\frac{\frac{g}{f}+\frac{1}{(1-f)}}{(g-1)^{2}}}=0$
(4) $0=\frac{\partial}{\partial
f}\sqrt{\frac{\frac{g}{f}+\frac{1}{(1-f)}}{(g-1)^{2}}}$ (5)

which, given $0<f<1$, is true for

$f=\frac{g-\sqrt{g}}{g-1}$ (6)

We can calculate the relative reduction in uncertainty compared to the 50/50
case (i.e. equal time spent on background and sample measurements) as:

$\frac{\frac{\Delta C_{bs}}{C_{bs}}\mid_{\mathrm{50/50}}-\frac{\Delta
C_{bs}}{C_{bs}}\mid_{\mathrm{optimal}}}{\frac{\Delta
C_{bs}}{C_{bs}}\mid_{\mathrm{50/50}}}=\frac{\sqrt{\frac{2g+2}{(g-1)^{2}}}-\sqrt{\frac{\frac{g^{2}-g}{g-\sqrt{g}}+\frac{g-1}{\sqrt{g}-1}}{(g-1)^{2}}}}{\sqrt{\frac{2g+2}{(g-1)^{2}}}}$
(7)

Figure 1: The optimal fraction of time $f$ spent measuring the sample as
opposed to measuring the background, as a function of the signal-to-noise
ratio $g$. Figure 2: The reduction in error that can be obtained by dividing
the time optimally between sample and background measurement, as a function of
the signal-to-noise ratio $g$.

##  III Conclusions

Figures 1 and 2 show the optimal division of time between sample and
background, and the reduction in uncertainty obtained through this
optimization, respectively. These figures clarify that the reduction in
uncertainty may be worth the trouble of a quick determination of the signal-
to-noise ratio, especially in areas where this ratio strongly deviates from
unity.

A quick scan of sample and background may be used to automatically select the
most optimal use of measurement time, in particular for scanning (small- and
wide-angle) diffractometers (including oddities such as Bonse-Hart cameras)
where the measurement time _per point_ can be freely tuned.

## References

  * Steinhart and Plestil (1993) M. Steinhart and J. Plestil, Journal of Applied Crystallography 26, 591 (1993).