Composite octet baryons in a relativistic mean field description of nuclear and neutron star matter

Kaito Noro 1ctad005@mail.u-tokai.ac.jp Graduate School of Science and Technology, Tokai University, 4-1-1 Kitakaname, Hiratsuka-shi, Kanagawa 259-1292, Japan Micro/Nano Technology Center, Tokai University, 4-1-1 Kitakaname, Hiratsuka-shi, Kanagawa 259-1292, Japan Wolfgang Bentz bentz@keyaki.cc.u-tokai.ac.jp Department of Physics, School of Science, Tokai University, 4-1-1 Kitakaname, Hiratsuka-shi, Kanagawa 259-1292, Japan Ian C. Cloët icloet@anl.gov Physics Division, Argonne National Laboratory, Argonne, Illinois 60439, USA Teruyuki Kitabayashi teruyuki@keyaki.cc.u-tokai.ac.jp Department of Physics, School of Science, Tokai University, 4-1-1 Kitakaname, Hiratsuka-shi, Kanagawa 259-1292, Japan

Abstract

We examine the properties of composite octet baryons in the nuclear medium and neutron star matter. The internal quark-diquark structure of the octet baryons and the equations of state of nuclear matter and neutron star matter in the mean field approximation are described by using the three-flavor Nambu–Jona-Lasinio (NJL) model as an effective quark theory of QCD. After introducing our model, we first discuss the properties of single baryons and their effective meson exchange interactions in symmetric nuclear matter by using concepts of Fermi liquid theory. Several model independent implications of this description are derived, and illustrated by numerical results obtained in our model. Second, we extend the model description to high baryon densities, and investigate the equation of state of neutron star matter and the resulting star masses. We find that the so called hyperon puzzle persists also for the case of composite hadrons. To get more information on this point, we also investigate the role of 6-fermi and 8-fermi interactions, in addition to the standard 4-fermi interactions. The strengths of those higher order fermi interactions is determined so as not to spoil the saturation properties of nuclear matter. Among them, an interaction characterized by a product of four quark current operators plays a special role to stabilize the stars over a large region of central baryon densities, although it has little effect on the maximum star masses. PhySH: Quark model; Asymmetric nuclear matter; Nuclear matter in neutron stars.

I INTRODUCTION

Systems of strongly interacting baryons are fascinating objects of current research, because their properties reflect the basic interaction between baryons which is intimately related to their quark substructure, and they connect microscopic nuclear systems to macroscopic astrophysical objects like supernovae and neutron stars. Besides the familiar building blocks of nuclear systems — protons and neutrons made of up ( $u$ ) and down ( $d$ ) quarks — baryons which carry strangeness are receiving much attention now because experimental and theoretical tools have become available to study their interactions and their role in nuclear and neutron star matter.

On the theoretical side, the baryon-baryon interactions have been extensively studied by using the meson exchange picture Haidenbauer and Meissner (2005); Rijken and Yamamoto (2006), effective field theories Bogner et al. (2010); Petschauer et al. (2020), and quantum Monte Carlo calculations Lonardoni et al. (2014). The parameters characterizing the two-body and possible three-body interactions are usually adjusted to scattering data, quark model predictions, or experimental data on nuclei and hypernuclei. Another line of approach, based on nonrelativistic constituent quark models, has been pursued vigorously Faessler et al. (1982); Oka et al. (1987); Fujiwara et al. (2007), mainly to understand the origin of the short-range repulsion on the basis of the Pauli principle on the quark level Oka (2023). More recent investigations are based on first principles derived from QCD Inoue et al. (2010); Hyodo and Niiyama (2021). These approaches, together, provide vital information to understand the properties of baryonic systems, in particular hypernuclei Gibson et al. (2010); Gal et al. (2016); Hiyama and Nakazawa (2018); Tamura (2022), and are useful tools to analyze new data on hyperon-nucleon scattering Miwa et al. (2022); Nanamura et al. (2022).

A test stone for theoretical models was provided by the observation of heavy neutron stars with about 2 solar masses Demorest et al. (2010); Antoniadis et al. (2013); Riley et al. (2021); Fonseca et al. (2021). Because the presence of hyperons usually leads to a softening of the equation of state of neutron star matter Glendenning (1997), many models were and are still unable to reproduce such heavy stars, and this problem is commonly called the “hyperon puzzle” Bombaci (2017). For extensive reviews on this subject and possible solutions, see for example Refs. Chatterjee and Vidaña (2016); Burgio et al. (2021). Most of the proposed solutions require additional repulsion between the baryons in the system, either via the exchange of vector mesons with particular forms of their couplings to baryons Weissenborn et al. (2012); Spinella and Weber (2019), pomeron exchange Yamamoto et al. (2013), or new kinds of three-body interactions Haidenbauer et al. (2017); Kohno (2018); Logoteta et al. (2019); Gerstung et al. (2020). Another possible solution Contrera et al. (2022) is based on the idea of a phase transition from nuclear matter to color superconducting quark matter Buballa (2005); Alford et al. (2008) at densities below or near the hyperon threshold.

As we mentioned at the beginning of this section, the properties of baryons and their interactions reflect their quark substructure, which changes in the nuclear medium. In order to study this aspect of the problem over a wide range of densities, relativistic quark models based on QCD are very useful tools. Two models of this kind, which have been used to describe nuclear phenomena in terms of quark degrees of freedom, are the quark-meson-coupling (QMC) model Guichon (1988), which is based on the MIT bag model Chodos et al. (1974), and the Nambu–Jona-Lasinio (NJL) model Nambu and Jona-Lasinio (1961a, b); Vogl and Weise (1991); Hatsuda and Kunihiro (1994), for which a full Faddeev approach Ishii et al. (1995) and a closely related but much simpler quark-diquark approach Bentz and Thomas (2001) to baryons have been developed. The degrees of freedom in the QMC model are quarks coupled to elementary mesons via Yukawa couplings, while the NJL model in its original form uses 4-fermi interactions between quarks to generate mesons as quark-antiquark bound states. Both models have been used extensively to explore the effects of medium modification on the quark level to nuclear observables Geesaman et al. (1995); Lu et al. (1998); Stone et al. (2016); Cloët et al. (2005); Cloet et al. (2006); Cloët et al. (2016).

The QMC model has also been applied to a wide range of hypernuclei Saito et al. (2007). Because meson exchange interactions usually tend to overbind the $\Lambda$ and $\Sigma$ baryons in nuclei Sammarruca (2008); Petschauer et al. (2020), in these earlier calculations a phenomenological repulsive interaction was introduced in order to reproduce the data. In a later version Guichon et al. (2008), the observation was made that the effect of spin-spin correlations between quarks, associated with the hyperfine interaction from gluon exchange, become enhanced in the nuclear medium if the $u$ and $d$ quark masses decrease as functions of the density but the strange ( $s$ ) quark mass remains constant. (For a simplified argument, see also Ref. Close (1979).) Because of the different spin-flavor structures of the $\Lambda$ and $\Sigma$ baryons, this leads to the expectation that the $\Sigma-\Lambda$ mass difference increases with the nuclear density. This kind of mechanism relies on the assumption of a constant $s$ quark mass, which is well satisfied in hypernuclei where the density of $s$ quarks is essentially zero, but may become less effective in neutron star matter as soon as a finite density of strange baryons appears.

In the present work, we will use the NJL model to describe the internal quark-diquark structure of the octet baryons, the equation of state of nuclear and neutron star matter in the mean field approximation, the corresponding in-medium effective meson exchange interactions between the baryons, and the resulting neutron star masses. The purposes of our work are as follows: First, we wish to explore the role of the quark-diquark substructure of baryons in the nuclear medium. For this purpose, we extend our previous work Carrillo-Serrano et al. (2016) on the properties of octet baryons in free space. Our model is well suited to examine the above expectation about the in-medium $\Sigma-\Lambda$ mass difference, because the spin-spin correlations in the scalar ( $0^{+}$ ) and axial vector ( $1^{+}$ ) diquark channels are built in from the outset. Second, in close connection to this, we wish to introduce ideas of the successful theory of Fermi liquids due to Landau Landau (1956, 1957, 1959) and Migdal Migdal (1967); Migdal et al. (1990), and its relativistic extensions Baym and Chin (1976), to hyperons in the nuclear medium. Because the power of the Fermi liquid theory to respect symmetries, conservation laws, and the renormalization group in many-fermion systems is well known Nozières (1964); Negele and Orland (1998); Shankar (1994), we find it desirable and timely to provide such a connection. Third, we wish to present a consistent formulation of isospin asymmetric baryonic systems on the background of the three independent Lorentz scalar and Lorentz vector mean fields, which are defined in Eq. (2) of the following section. Finally, we wish to investigate the status of the hyperon puzzle in the NJL model, and investigate the roles of 6-fermi ’t Hooft (1976) and 8-fermi Osipov et al. (2007) interactions on the equation of state and star masses in the mean field approximation. In order to achieve these aims as clearly as possible, we will make no attempt to reproduce any empirical data related to octet baryons, their mutual interactions, or properties of neutron stars. Rather than this, we wish to explain problems which arise from chiral symmetry restrictions on the form of the interaction Lagrangian, which were not encountered in our previous work on the flavor SU(2) case Tanimoto et al. (2020).

The outline of the paper is as follows: Sec. II discusses our effective quark model for octet baryons and baryonic matter; Sec. III discusses the properties of baryons and their effective meson exchange interactions in symmetric nuclear matter using concepts of Fermi liquid theory; Sec. IV presents our results for neutron star matter and the resulting star masses; Sec. V discusses the roles of 6-fermi and 8-fermi interactions; and Sec. VI gives a summary of our results.

II MODEL FOR BARYONS AND BARYONIC MATTER

The three-flavor NJL Lagrangian with 4-fermi interactions in the $\bar{q}q$ channels relevant for this study reads Vogl and Weise (1991); Hatsuda and Kunihiro (1994):

	$\displaystyle\mathcal{L}$	$\displaystyle=\bar{q}\left(i\not{\partial}-\hat{m}\right)q+G_{\pi}\left[\big{(}\bar{q}\lambda_{a}q\big{)}^{2}-\big{(}\bar{q}\lambda_{a}\gamma_{5}q\big{)}^{2}\right]\allowdisplaybreaks$
		$\displaystyle-G_{v}\left[\big{(}\bar{q}\lambda_{a}\gamma^{\mu}q\big{)}^{2}+\big{(}\bar{q}\lambda_{a}\gamma^{\mu}\gamma_{5}q\big{)}^{2}\right]\,,$		(1)

where $q=(q_{1},q_{2},q_{3})$ with $1\equiv u$ , $2\equiv d$ , $3\equiv s$ is the quark field, $\hat{m}$ the current quark mass matrix with diagonal elements $(m_{u},m_{d},m_{s})$ , and $\lambda_{a}$ ( $a=0,1,2,\dots 8$ ) are the Gell-Mann flavor matrices plus $\lambda_{0}=\sqrt{\frac{2}{3}}\mathbf{1}$ . The 4-fermi coupling constants in the scalar–pseudoscalar and the vector–axial-vector channels are denoted by $G_{\pi}$ and $G_{v}$ , respectively. The Lagrangian (1) has the $SU(3)_{L}\otimes SU(3)_{R}\otimes U(1)_{V}\otimes U(1)_{A}$ symmetry of QCD, which contains the familiar flavor $SU(3)$ as a subgroup. The explicit breaking of the $U(1)_{A}$ symmetry, which is known as the axial anomaly in QCD, can be realized in the NJL model by the 6-fermi (determinant) interaction ’t Hooft (1976), which will be investigated together with possible 8-fermi interactions in Sec. V. It is important to note that in this work we will follow the successful path established by various low energy theorems and octet mass formulas, that current quark masses are the only sources of explicit breaking of the flavor and the chiral symmetries, and all other symmetry breakings are dynamical. As we will see, this leads to very strong, sometimes unwelcome, restrictions on the model parameters in the mean field approximation.

In order to construct the octet baryons as quark-diquark bound states, we also need the interaction Lagrangian in the $qq$ channels with the same symmetries, which is specified in App. A. Our model description of the octet baryons is a straight forward extension of the quark-diquark model based on the Faddeev framework, as described in Refs. Carrillo-Serrano et al. (2014, 2016), to the case where the isospin symmetry is broken, like in neutron star matter. In the vacuum isospin symmetry is assumed to be intact, i.e., we use $m_{u}=m_{d}\equiv m$ throughout this work.

II.1 Mean field approximation

In order to construct the equation of state of nuclear matter and neutron star matter in the mean field approximation, we will take into account three scalar fields $\sigma_{\alpha}$ and three 4-vector fields $\omega_{\alpha}^{\mu}$ , where $\alpha=u,d,s$ . We use the following definitions:

\displaystyle\sigma_{\alpha}

\displaystyle=4G_{\pi}\langle\bar{q}_{\alpha}\,q_{\alpha}\rangle\,,

\displaystyle\omega_{\alpha}^{\mu}

\displaystyle=4G_{v}\langle\bar{q}_{\alpha}\gamma^{\mu}\,q_{\alpha}\rangle\,,

(2)

where $\langle\dots\rangle$ denotes the expectation value in the ground state of the medium under consideration (vacuum, nuclear matter, or neutron star matter). The presence of the scalar fields leads to spontaneous breaking of the chiral symmetry, and gives rise to the effective quark masses

\displaystyle M_{\alpha}=m_{\alpha}-\sigma_{\alpha}\,,

(3)

which must be treated independently if the isospin symmetry is broken in the medium. The presence of the vector fields leads to shifts in the 4-momenta of the particles in the system. As a result, the energy of a baryon with flavor $b$ and 3-momentum $\bm{k}$ is obtained from the pole of the quark-diquark equation in the variable $k_{0}$ as¹¹1Here and in the following, a summation over multiple flavor indices ( $\alpha,\beta,\dots$ for quarks, $b,b^{\prime}$ for octet baryons, $\tau$ for the special case of nucleons, and $i$ for baryons and leptons) in a product, including squares like $\omega_{\alpha}^{2}$ , is implied if those indices appear only on one side of an equation. (As usual, the same convention is used for the Lorentz indices $\mu,\nu,\dots$ .) The Fermi momentum of particle $i$ will be denoted as $p_{i}$ .

\displaystyle\varepsilon_{b}(k)=\sqrt{\bm{k}_{b}^{2}+M_{b}^{2}}+n_{\alpha/b}\,\omega^{0}_{\alpha}\equiv E_{b}(k_{b})+n_{\alpha/b}\,\omega^{0}_{\alpha}\,,

(4)

where $n_{\alpha/b}$ is the number of quarks with flavor $\alpha$ in the baryon $b$ , and $\bm{k}_{b}=\bm{k}-n_{\alpha/b}\,\bm{\omega}_{\alpha}$ . The effective mass of the baryon, $M_{b}$ , is a function of the effective quark masses $M_{u},M_{d},M_{s}$ , as described in App. A.

The mean field approximation is implemented into the Lagrangian (1) in the standard way by decomposing the various quark bilinears into classical (c-number) parts and quantum (normal ordered) parts. We will assume that the only non-vanishing classical parts are the mean fields given in Eq. (2). The normal ordered parts, together with the $qq$ interaction parts given in App. A, are used to calculate bound state masses of pseudoscalar mesons and octet baryons, as well as the pion decay constant.

The quantity of central interest in our work is the energy density ( ${\cal E}$ ) of baryonic matter in the mean field approximation. The basic physical picture can be visualized by composite baryons moving in scalar and vector mean fields on the background of the constituent quark vacuum Bentz and Thomas (2001). Except for the vacuum contributions, this is similar in spirit to the QMC model Guichon (1988); Saito and Thomas (1994), although the mesons in our approach are composite objects. The term which describes the Fermi motion of the baryons is given by (note our summation convention for multiple flavor indices)

	$\displaystyle 2\int\frac{{\rm d}^{3}k}{(2\pi)^{3}}\,\varepsilon_{b}(k)\,n_{b}(k)$
	$\displaystyle=2\int\frac{{\rm d}^{3}k}{(2\pi)^{3}}E_{b}(k)\,n_{b}(k)+\rho_{\alpha}\,\omega^{0}_{\alpha}\equiv{\cal E}_{B}+\rho_{\alpha}\,\omega^{0}_{\alpha}\,,$		(5)

where $n_{b}(k)$ is the Fermi distribution function of baryon $b$ , and we defined the quark number densities $\rho_{\alpha}$ in terms of the baryon number densities $\rho_{b}$ by $\rho_{\alpha}=n_{\alpha/b}\rho_{b}$ . For the case of neutron star matter we also include the contributions from the Fermi gas of leptons ( $\ell=e^{-},\mu^{-}$ ) in chemical equilibrium with the baryons. The total energy density in the mean field approximation is then expressed as

\displaystyle{\cal E}

\displaystyle={\cal E}_{\rm vac}-\frac{\omega_{\alpha}^{2}}{8G_{v}}+\rho_{\alpha}\omega^{0}_{\alpha}+{\cal E}_{B}+{\cal E}_{\ell}\,.

(6)

Here the unregularized form of the vacuum (Mexican hat shaped) contribution is

\displaystyle{\cal E}_{\rm vac}=6i\int\frac{{\rm d}^{4}k}{(2\pi)^{4}}\,\ln\frac{k^{2}-M_{\alpha}^{2}}{k^{2}-M_{\alpha 0}^{2}}+\frac{\sigma_{\alpha}^{2}-\sigma_{\alpha 0}^{2}}{8G_{\pi}}\,,

(7)

where a sum over the quark flavors $\alpha$ is implied, and the sub-index $0$ refers to the vacuum with zero baryon density.

The scalar and vector fields are determined for given baryon density $\rho_{B}$ by the conditions

\displaystyle\partial{\cal E}/(\partial\sigma_{\alpha})=\partial{\cal E}/(\partial\omega^{\mu}_{\alpha})=0\,.

(8)

For the scalar fields, the minimizations (8) have to be done numerically. It is, however, easy to confirm that they are equivalent to the relation

\displaystyle\sigma_{\alpha}=4\,G_{\pi}\frac{\partial{\cal E}}{\partial m_{\alpha}}=4\,G_{\pi}\,\frac{\partial{\cal E}}{\partial M_{\alpha}}\,,

(9)

where the first equality is the general Feynman-Hellman theorem, while the second equality holds if the energy density is expressed in such a way that the constituent quark masses $M_{\alpha}$ always appear together with the current quark masses $m_{\alpha}$ , i.e., in the first term of the vacuum energy (7) and in the term ${\cal E}_{B}$ of (5). For the vector fields, Eq. (8) leads to

\displaystyle\omega_{\alpha}^{\mu}=4\,G_{v}\,j_{\alpha}^{\mu}=4\,G_{v}\,n_{\alpha/b}\,j_{b}^{\mu}\,,

(10)

where $j_{\alpha}^{\mu}=\left(\rho_{\alpha},\bm{j}_{\alpha}\right)$ is the contribution to the baryon current carried by the quark of flavor $\alpha$ , and $j_{b}^{\mu}$ is the corresponding quantity for the baryon $b$ . Eq. (10) is in accordance with the definition given in Eq. (2).

For neutron star matter, the minimization w.r.t. the scalar fields — or equivalently the solution to Eq. (9) — has to be done under the requirements of chemical equilibrium and charge neutrality Glendenning (1997)

\displaystyle\mu_{b}-\mu_{n}+q_{b}\,\mu_{e}=\mu_{\mu}-\mu_{e}=q_{i}\,\rho_{i}=0\,,

(11)

where the chemical potentials for baryons and leptons are given by $\mu_{b}=\varepsilon_{b}(k=p_{b})$ and $\mu_{\ell}=\sqrt{p_{\ell}^{2}+m_{\ell}^{2}}$ . The Fermi momenta $p_{i}$ of baryons ( $i=b$ ) and leptons ( $i=\ell$ ) are related to their number densities $\rho_{i}$ by $\rho_{i}=\frac{p_{i}^{3}}{3\pi^{2}}$ . In Eq. (11), $q_{i}$ are the electric charges of baryons and leptons. In the general case, for given baryon density, the nine independent relations in Eq. (11) determine the densities of 10 particles in the system (8 baryons and 2 leptons). The pressure of the system can then be obtained as a function of baryon density from the relation $P=\rho_{i}\,\mu_{i}\,-{\cal E}$ .

II.2 Effective baryon-baryon interaction

For the purpose of discussions, it will be useful to know the form of the effective baryon-baryon interaction which underlies the mean field approximation described above. For this purpose, we follow the ideas of the Fermi liquid theory Nozières (1964); Migdal (1967) and its relativistic extensions Baym and Chin (1976), and define the spin averaged effective baryon-baryon interaction $F_{bb^{\prime}}(p,p^{\prime})$ by the variation of the energy of one of the baryons, $\varepsilon_{b}(p)$ , w.r.t. the distribution function of the other baryon, $n_{b^{\prime}}(p^{\prime})$ . We wish to express this interaction as a generalized meson-exchange potential. Because our baryon energies in Eq. (4) do no depend explicitly on the distribution functions, we have

\displaystyle F_{bb^{\prime}}=\frac{\delta\varepsilon_{b}}{\delta n_{b^{\prime}}}=\frac{\partial\varepsilon_{b}}{\partial\sigma_{\alpha}}\frac{\delta\sigma_{\alpha}}{\delta n_{b^{\prime}}}+\frac{\partial\varepsilon_{b}}{\partial\omega^{\mu}_{\alpha}}\frac{\delta\omega^{\mu}_{\alpha}}{\delta n_{b^{\prime}}}\,,

(12)

where we omitted the dependence on the momenta $p$ and $p^{\prime}$ to simplify the notations. Because the conditions given in Eq. (8) hold for any fixed set of distribution functions, we can make use of the relations

	$\displaystyle\frac{\delta}{\delta n_{b^{\prime}}}\left(\frac{\partial{\cal E}}{\partial\sigma_{\alpha}}\right)$	$\displaystyle=0=\frac{\partial{\varepsilon_{b^{\prime}}}}{\partial\sigma_{\alpha}}+\frac{\partial^{2}{\cal E}}{\partial\sigma_{\alpha}\partial\sigma_{\beta}}\,\frac{\delta\sigma_{\beta}}{\delta n_{b^{\prime}}}\,,$
	$\displaystyle\frac{\delta}{\delta n_{b^{\prime}}}\left(\frac{\partial{\cal E}}{\partial\omega^{\mu}_{\alpha}}\right)$	$\displaystyle=0=\frac{\partial{\varepsilon_{b^{\prime}}}}{\partial\omega^{\mu}_{\alpha}}+\frac{\partial^{2}{\cal E}}{\partial\omega^{\mu}_{\alpha}\partial\omega^{\nu}_{\beta}}\,\frac{\delta\omega^{\nu}_{\beta}}{\delta n_{b^{\prime}}}\,,$		(13)

where the second equalities hold in our model when the whole system is at rest, in which case there are no mixings between scalar and vector mean fields. Using (13) in (12) we obtain

	$\displaystyle F_{bb^{\prime}}$	$\displaystyle=-\frac{M_{b}}{E_{b}}\,\frac{\partial M_{b}}{\partial\sigma_{\alpha}}\left(S^{-1}\right)_{\alpha\beta}\frac{\partial M_{b^{\prime}}}{\partial\sigma_{\beta}}\,\frac{M_{b^{\prime}}}{E_{b^{\prime}}}-n_{\alpha/b}\left(V^{-1}\right)^{00}_{\alpha\beta}\,n_{\beta/b^{\prime}}$
		$\displaystyle-\frac{p_{i}}{E_{b}}\,n_{\alpha/b}\left(V^{-1}\right)^{ij}_{\alpha\beta}n_{\beta/b}\,\frac{p^{\prime}_{j}}{E_{b^{\prime}}}\,.$		(14)

Here $E_{b}\equiv E_{b}(p)$ , $E_{b^{\prime}}\equiv E_{b^{\prime}}(p^{\prime})$ , and we defined the $3\times 3$ flavor matrices $S$ and $V$ by

\displaystyle S_{\alpha\beta}

\displaystyle\equiv\frac{\partial^{2}{\cal E}}{\partial\sigma_{\alpha}\partial\sigma_{\beta}}\,,

\displaystyle V^{\mu\nu}_{\alpha\beta}

\displaystyle\equiv\frac{\partial^{2}{\cal E}}{\partial\omega_{\alpha\mu}\partial\omega_{\beta\nu}}\,.

(15)

We illustrate the effective interaction of Eq. (14) by Fig. 1, where the solid lines express the baryons, the dashed line expresses the generalized propagators of neutral scalar mesons ( $S^{-1}$ ) and vector meson ( $V^{-1}$ ) for zero momenta, and the vertices stand for the factors to the left and the right of the meson propagators in Eq. (14).

In isospin asymmetric baryonic matter, like neutron star matter, the $\bar{u}u$ , $\bar{d}d$ and $\bar{s}s$ components of the exchanged mesons are mixed by the baryon loop term ${\cal E}_{B}$ . To disentangle them, one could make an orthogonal transformation to diagonalize $S$ and $V$ at fixed baryon density, and express the couplings of each exchanged flavor to the baryon by a linear combination of vertices. In the present work we will not carry out such an analysis for the case of neutron star matter. In the case of isospin symmetric nuclear matter, on the other hand, the matrices $S$ and $V$ become diagonal automatically by taking isoscalar and isovector combinations of the interacting baryons in the particle-hole channel ( $t$ -channel), and we will show the explicit forms in the next section.

Refer to caption — Figure 1: Graphical representation of the effective baryon-baryon interaction (14) as a meson exchange potential. For explanation of symbols, see the text.

III BARYONS IN SYMMETRIC NUCLEAR MATTER

In this section we wish to discuss our results for the properties of baryons and their mutual effective meson exchange interactions in isospin symmetric nuclear matter. In this case, the mean fields (2) with $\alpha=u$ and $\alpha=d$ are the same because of the isospin symmetry, and $\omega_{s}^{\mu}=0$ because the density of strange quarks is zero. Because the $s$ -quark mass enters only in the vacuum energy (7), the minimization condition $\partial{\cal E}/\partial\sigma_{s}=0$ gives $\sigma_{s}=\sigma_{s0}$ and therefore $M_{s}=M_{s0}$ . (Note that this holds only in the present case of 4-fermi interactions. The 6-fermi and 8-fermi interactions considered in Sec. V lead to a slight density dependence of $M_{s}$ even in symmetric nuclear matter.) The energy density of the system is given by Eq. (6) without the leptonic term ${\cal E}_{\ell}$ .

III.1 Effective meson exchange interaction

In the case of isospin symmetric nuclear matter, the flavor matrices of Eq. (15), which characterize the effective interaction (14), become diagonal automatically by taking appropriate combinations, e.g., for the case of the baryon-nucleon interaction we define

\displaystyle f_{bN}

\displaystyle\equiv\frac{1}{2}\left(F_{bp}+F_{bn}\right),

\displaystyle f^{\prime}_{bN}

\displaystyle\equiv\frac{1}{2}\left(F_{bp}-F_{bn}\right).

(16)

Within the isospin multiplet to which the baryon $b$ belongs, $f_{bN}$ is an isoscalar and the same for all members of the multiplet, while $f^{\prime}_{bN}$ is an isovector proportional to the isospin 3-component of the baryon $b$ . We find the explicit forms

	$\displaystyle f_{bN}(p,p^{\prime})=-\frac{M_{b}}{E_{b}}\frac{M_{N}}{E_{N}}\left(\frac{\partial M_{b}}{\partial M}\right)\left(\frac{\partial M_{N}}{\partial M}\right)$
	$\displaystyle\times\left(\frac{1}{2G_{\pi}}+2g(M)+\rho_{N}^{(s)}\frac{\partial^{2}M_{N}}{\partial M^{2}}+\phi_{N}\left(\frac{\partial M_{N}}{\partial M}\right)^{2}\right)^{-1}$
	$\displaystyle+\frac{6\left(1+\frac{y_{b}}{2}\right)}{1/(2G_{v})}-6\left(1+\frac{y_{b}}{2}\right)\frac{\bm{p}\cdot\bm{p}^{\prime}}{E_{b}\,E_{N}}\,\left(\frac{1}{2G_{v}}+\frac{9\rho_{B}}{E_{N}}\right)^{-1},$		(17)

	$\displaystyle f^{\prime}_{bN}(p,p^{\prime})=-\frac{M_{b}}{E_{b}}\frac{M_{N}}{E_{N}}\left(\frac{\partial M_{b}}{\partial\Delta M}\right)\left(\frac{\partial M_{p}}{\partial\Delta M}\right)$
	$\displaystyle\times\left(\frac{1}{2G_{\pi}}+2g(M)+\rho_{N}^{(s)}\frac{\partial^{2}M_{p}}{\partial(\Delta M)^{2}}+\phi_{N}\left(\frac{\partial M_{p}}{\partial\Delta M}\right)^{2}\right)^{-1}$
	$\displaystyle+\frac{2\,t_{b}}{1/(2G_{v})}-2t_{b}\,\frac{\bm{p}\cdot\bm{p}^{\prime}}{E_{b}\,E_{N}}\,\left(\frac{1}{2G_{v}}+\frac{\rho_{B}}{E_{N}}\right)^{-1}.$		(18)

Here $M$ is the effective mass of $u,d$ quarks, $M_{N}$ the nucleon mass in symmetric nuclear matter, and the derivatives w.r.t. $\Delta M\equiv M_{u}-M_{d}$ should be evaluated at $M_{u}=M_{d}$ . The unregularized form of the function $g(M)$ is

\displaystyle g(M)=-12i\int\frac{{\rm d}^{4}k}{(2\pi)^{4}}\frac{k^{2}+M^{2}}{\left(k^{2}-M^{2}\right)^{2}}\,,

(19)

and the scalar density of the nucleon ( $\rho_{N}^{(s)}$ ) and the function $\phi_{N}$ are defined by

	$\displaystyle\rho_{N}^{(s)}$	$\displaystyle=2\int\frac{{\rm d}^{3}k}{(2\pi)^{3}}\,n_{N}(k)\,\frac{M_{N}}{E_{N}(k)}\,,$		(20)
	$\displaystyle\phi_{N}$	$\displaystyle=2\int\frac{{\rm d}^{3}k}{(2\pi)^{3}}\,n_{N}(k)\,\frac{k^{2}}{E_{N}(k)^{3}}\,.$		(21)

In Eqs. (17) and (18), $E_{b}\equiv E_{b}(p)$ , $E_{N}\equiv E_{N}(p^{\prime})$ , and we used the relations $n_{u/b}+n_{d/b}=2\left(1+\frac{y_{b}}{2}\right)$ and $n_{u/b}-n_{d/b}=2t_{b}$ where $y_{b}$ and $t_{b}$ are the hypercharge and the isospin 3-component of the baryon $b$ .

The interpretation of (17) for the case $b=N$ in terms of the meson exchange processes of Fig. 1 has been discussed in detail in Ref. Bentz and Thomas (2001), and the generalization is almost self evident:²²2To make the connection to the baryon-meson coupling constants and meson masses, the numerator functions and denominator functions (i.e., those parts which involve $\frac{1}{2G_{\pi}}$ or $\frac{1}{2G_{v}}$ ) must be multiplied by the squares of the relevant quark-meson coupling constants, i.e., by $g^{(q)2}_{\sigma}=g^{(q)2}_{\delta}$ for the first 2 lines of (17) and (18), and by $g^{(q)2}_{\omega}=g^{(q)2}_{\rho}$ for the third lines. See App. B for details. The first two lines (third line) in (17) correspond to $\sigma$ ( $\omega$ ) meson exchange, while the first two lines (third line) in (18) correspond to neutral $\delta$ ( $\rho$ ) meson exchange. The coupling constants of a baryon $b$ to the $\sigma$ ( $\delta$ ) meson are proportional to the derivative of $M_{b}$ w.r.t. $M$ ( $\Delta M$ ), while the couplings to the vector mesons ( $\omega$ and $\rho$ ) are determined by the isoscalar and isovector combination of the quark numbers in the baryon. The function $g(M)$ in the denominators of (17) and (18) is the one quark-loop self energy of the scalar meson in the vacuum, the terms involving the scalar density are the Fermi averages over effective $\sigma\sigma NN$ ( $\delta\delta NN$ ) contact interactions which are induced by the scalar-isoscalar polarizability $\frac{\partial^{2}M_{N}}{\partial M^{2}}$ (the scalar-isovector polarizability $\frac{\partial^{2}M_{p}}{\partial(\Delta M)^{2}}$ ) of the nucleon Birse (1995); Saito and Thomas (1995); Wallace et al. (1995), and the terms proportional to $\phi_{N}$ are the Fermi averages over the “Z-graph” contributions, which also appear in hadronic theories Brown and Weise (1987); Wallace (1998). In our numerical calculations, discussed in Sec. III.3, we find that the numerators of the $\sigma$ and $\delta$ exchange parts in (17) and (18) substantially decrease, while their denominators slightly increase as the baryon density increases. The increase of the denominators is partially related to the fact that the scalar-isoscalar and scalar-isovector polarizabilities of nucleons — both being positive — increase as the baryon density increases. The density dependence of the coupling constants and mesons self energies in (17) and (18) then suppresses the attractive effects of scalar meson exchanges, and at higher densities the vector meson exchanges become dominant. The terms $\propto\bm{p}\cdot\bm{p}^{\prime}$ correspond to the contributions from the exchange of $\bm{\omega}$ and neutral $\bm{\rho}$ mesons, and their self energies arise only from the Fermi averages over the corresponding Z-graphs. More details will be given in Sec. III.3.3 and App. B.

III.2 Physical implications of the interaction

In order to explain some physical implications of the effective meson exchange interactions $f_{bN}$ and $f^{\prime}_{bN}$ of Eqs. (17) and (18), we extend a few basic points of Fermi liquid theory to octet baryons in the nuclear medium. In the following discussions, we will use the following notations:³³3We remind that “isovector” ( $T=1$ ) and “isoscalar” ( $T=0$ ) refers to the particle-hole channel (t-channel) of the interacting baryons, not to the isospin of the two incoming baryons.

•

$B=N,\Sigma,\Lambda,\Xi$ stands for the isospin multiplets (including the isospin singlet $\Lambda$ ), while $b$ continues to stand for a member of the baryon octet;
•

$f_{BN}$ denotes the 4 independent isoscalar baryon-nucleon interactions, defined by Eq. (17) with $b$ a member of $B$ ;
•

$f^{\prime}_{BN}$ denotes the 3 independent isovector baryon-nucleon interactions ( $B=N,\Sigma,\Xi$ ), defined by Eq. (18) with $b$ a member of $B$ with the largest value of the isospin 3-component $t_{b}$ (i.e., $p$ , $\Sigma^{+}$ , and $\Xi^{0}$ ).

We also separate the terms $\propto\bm{p}\cdot\bm{p}^{\prime}$ , which involve the transfer of one unit of orbital angular momentum ( $\ell=1$ ) between the baryons, from the other terms which involve no angular momentum transfer ( $\ell=0$ ):

\displaystyle f_{bN}

\displaystyle=f_{0,bN}+\cos\theta\,f_{1,bN}\,,

\displaystyle f^{\prime}_{bN}

\displaystyle=f^{\prime}_{0,bN}+\cos\theta\,f^{\prime}_{1,bN}\,,

(22)

and similarly for $f_{BN}$ , $f^{\prime}_{BN}$ , $F_{bp}$ , and $F_{bn}$ , where $\theta$ is the angle between $\bm{p}$ and $\bm{p}^{\prime}$ . The parameters $f_{\ell,NN}$ and $f^{\prime}_{\ell,NN}$ defined in this way agree with the familiar Landau-Migdal parameters, usually denoted by $f_{\ell}$ and $f^{\prime}_{\ell}$ .

III.2.1 Nucleon density variations

The $\ell=0$ baryon-nucleon interactions $f_{0,bN}$ and $f^{\prime}_{0,bN}$ express the change of the baryon energies, $\varepsilon_{b}$ , caused by variations of the Fermi momenta of the background nucleons. If we denote the corresponding variations of nucleon densities by $\delta\rho_{\tau}$ ( $\tau=p,n$ ), the change of the distribution functions to first order is given by $\delta n_{\tau}(k)=\frac{\pi^{2}}{p^{2}_{\tau}}\,(\delta\rho_{\tau})\,\delta(k-p_{\tau})$ . Then, according to the general definition given by the first equality in Eq. (12), the energy of a baryon $b$ in nuclear matter changes by an amount

\displaystyle\delta\varepsilon_{b}(k)

\displaystyle=2\int\frac{{\rm d}^{3}p}{(2\pi)^{3}}\,F_{b\tau}(k,p)\,\delta n_{\tau}(p)=\,F_{0,b\tau}\,\delta\rho_{\tau}\,.

(23)

Separating the isoscalar from the isovector contributions then gives the general relations

\displaystyle\frac{\partial\varepsilon_{b}(k)}{\partial\rho_{B}}

\displaystyle=f_{0,bN}(k,p_{N})\,,

\displaystyle\frac{\partial\varepsilon_{b}(k)}{\partial\rho_{(3)}}

\displaystyle=f^{\prime}_{0,bN}(k,p_{N})\,,

(24)

where $\rho_{(3)}=\rho_{p}-\rho_{n}$ , and the limit of isospin symmetric nuclear matter ( $\rho_{(3)}\rightarrow 0$ ) is understood. For the case where $b$ is a nucleon, the parameters $f_{0,NN}$ and $f^{\prime}_{0,NN}$ for $k=p_{N}$ are related to the incompressibility ( $K$ ) and the symmetry energy ( $a_{s}$ ) as follows:

	$\displaystyle K$	$\displaystyle=9\rho_{B}\left(\frac{\pi^{2}}{2p_{N}E_{N}}+f_{0,NN}\right),$		(25)
	$\displaystyle a_{s}$	$\displaystyle=\frac{\rho_{B}}{2}\left(\frac{\pi^{2}}{2p_{N}E_{N}}+f^{\prime}_{0,NN}\right).$		(26)

In our model the baryon energy is given by Eq. (4), and by using Eq. (10) in nuclear matter at rest, we have

\displaystyle\varepsilon_{b}(k)=E_{b}(k)+12G_{v}\rho_{B}\,\left(1+\frac{y_{b}}{2}\right)+4G_{v}\,\rho_{(3)}\,t_{b}\,.

(27)

Using this in Eq. (24), we see that in our model the $\ell=0$ baryon-nucleon interaction reflects the density dependence of the baryon effective masses:

	$\displaystyle f_{0,bN}(k,p_{N})$	$\displaystyle=\frac{M_{b}}{E_{b}(k)}\,\frac{\partial M_{b}}{\partial\rho_{B}}+12G_{v}\left(1+\frac{y_{b}}{2}\right)\,,$		(28)
	$\displaystyle f^{\prime}_{0,bN}(k,p_{N})$	$\displaystyle=\frac{M_{b}}{E_{b}(k)}\,\frac{\partial M_{b}}{\partial\rho_{(3)}}+4\,G_{v}\,t_{b}\,.$		(29)

For the case where $b$ is a nucleon, the two terms in $f_{0,NN}$ , when multiplied by $9\rho_{B}$ , give the contributions of $\sigma$ meson and $\omega$ meson exchange to the incompressibility. Similarly, the two terms in $f^{\prime}_{0,NN}$ , when multiplied by $\frac{\rho_{B}}{2}$ , give the contributions of $\delta$ meson and $\rho$ meson exchange to the symmetry energy.

III.2.2 Lorentz invariance

There are two basic requirements from Lorentz invariance in the present context: First, the distribution function of the nucleons is Lorentz invariant: $n^{\prime}_{\tau}(k^{\prime})=n_{\tau}(k)$ , where we use a prime to denote a system which moves with velocity $\bm{u}$ relative to the reference system which we assume to be at rest, and $k^{\prime}=\Lambda_{\bm{u}}\,k$ where $\Lambda_{\bm{u}}$ is the Lorentz matrix. A Lorentz transformation then leads to a variation of the distribution function for fixed momentum according to Baym and Chin (1976) $\delta n_{\tau}(k)\equiv n_{\tau}^{\prime}(k)-n_{\tau}(k)=-\varepsilon_{\tau}(k)\bm{u}\cdot\hat{\bm{k}}\delta(k-p_{\tau})$ to first order in $\bm{u}$ . Second, the change of the energy of a baryon in symmetric nuclear matter, induced by this density variation,

	$\displaystyle\delta\varepsilon_{b}(k)$	$\displaystyle=2\int\frac{{\rm d}^{3}p}{(2\pi)^{3}}\,F_{b\tau}(k,p)\,\delta n_{\tau}(p)$
		$\displaystyle=-\bm{u}\cdot\hat{\bm{k}}\,\frac{2p_{N}^{2}\varepsilon_{N}}{3\pi^{2}}\,f_{1,bN}(k,p_{N})\,,$		(30)

must be equivalent to a Lorentz transformation applied directly to the baryon energy, $\delta\varepsilon_{b}(k)\equiv\varepsilon_{b}^{\prime}(k)-\varepsilon_{b}(k)=-\bm{u}\cdot\bm{k}+\varepsilon_{b}(k)\,\bm{u}\cdot\bm{v}_{b}(k)$ to first order in $\bm{u}$ , where $\bm{v}_{b}(k)=\bm{\nabla}_{k}\varepsilon_{b}(k)$ is the velocity of the baryon. This requirement leads to the relation

\displaystyle\frac{k}{\varepsilon_{b}(k)}=v_{b}(k)+\frac{\varepsilon_{N}}{\varepsilon_{b}(k)}\frac{2p_{N}^{2}}{3\pi^{2}}\,f_{1,bN}(k,p_{N})\,.

(31)

In Eqs. (30) and (31), and in all following relations, $\varepsilon_{N}\equiv\varepsilon_{N}(p_{N})$ is the Fermi energy of the nucleon in symmetric nuclear matter, while the momentum $k$ is arbitrary. For the case where $b$ is a nucleon, Eq. (31) agrees with the relativistic form of the Landau effective mass relation for variable momentum $k$ Baym and Chin (1976). The velocity of the baryon is usually expressed in terms of the Landau effective mass ( $M_{b}^{*}(k)$ ) by $v_{b}(k)\equiv\frac{k}{M_{b}^{*}(k)}$ . By taking the limits $k\rightarrow 0$ on both sides of Eq. (31), we then obtain a simple relation of the form

\displaystyle\frac{1}{\varepsilon_{b}(0)}=\frac{1}{M_{b}^{*}(0)}+\frac{\varepsilon_{N}}{\varepsilon_{b}(0)}\,\rho_{B}\,\hat{f}_{1,bN}(0,p_{N})\,,

(32)

where we defined $\hat{f}_{1,bN}(p,p^{\prime})$ such that the “full” $\ell=1$ interaction (like for example the last term in Eq. (17)) is expressed in the form $\left(\bm{p}\cdot\bm{p}^{\prime}\right)\,\hat{f}_{1,bN}(p,p^{\prime})$ .

It is easy to check that our model satisfies the requirement (31): The energy of a baryon with momentum $k$ and the Fermi energy $\varepsilon_{N}$ of a nucleon in symmetric nuclear matter are obtained from (27) by setting $\rho_{(3)}=0$ :

	$\displaystyle\varepsilon_{b}(k)$	$\displaystyle=E_{b}(k)+12\,G_{v}\rho_{B}\left(1+\frac{y_{b}}{2}\right)\,,$
		$\displaystyle\varepsilon_{N}=E_{N}+18\,G_{v}\rho_{B}\,,$		(33)

while $f_{1,bN}$ , which corresponds to $\bm{\omega}$ exchange, is given by the last term in Eq. (17) without the factor $\cos\theta$ , see Eq. (22). It is then clear that the general relation (31) is valid in our model.

III.2.3 Currents carried by baryons

The Lorentz invariance requirement of Eq. (31) is related to the isoscalar $\ell=1$ Fermi liquid parameter $f_{1,bN}$ . To give an example where also the isovector part enters, let us consider the currents carried by a baryon $b$ moving with momentum $k$ in nuclear matter, for the case where no momentum is transferred by the external field to the baryon. From gauge invariance and the integral equations for the vertex functions, the Fermi liquid theory leads to the following result Bentz et al. (1985) (see also, for example, Eqs. (1)-(33) of Ref. Nozières (1964) or Eq. (2.16) of Ref. Migdal (1967)):

\displaystyle\bm{j}_{b}^{(X)}(k)

\displaystyle=\bm{v}_{b}(k)\,Q_{b}^{(X)}+2\,p_{\tau}^{2}\,Q^{(X)}_{\tau}\int\frac{{\rm d}\Omega_{p}}{(2\pi)^{3}}\,\hat{\bm{p}}\,F_{b\tau}(k,p)\,,

(34)

where $X$ characterizes the type of current, e.g., $X=B$ for the baryon current, $X=I$ for the isospin current, and $X=E$ for the electric current, and $Q_{b}^{(X)}$ are the corresponding bare charges of the baryon $b$ , i.e., $Q_{b}^{(B)}=1$ , $Q_{b}^{(I)}=t_{b}$ , and $Q_{b}^{(E)}=q_{b}$ . The second term in Eq. (34) is the backflow due to the nuclear medium.

The magnitude of the baryon current (case $X=B$ in Eq. (34)) can be expressed in a model independent way by using the Lorentz invariance relation (31):

\displaystyle j_{b}^{(B)}(k)

\displaystyle=\frac{k}{\varepsilon_{N}}-v_{b}(k)\left(\frac{\varepsilon_{b}(k)}{\varepsilon_{N}}-1\right)\,.

(35)

For the case of a nucleon at the Fermi surface, Eq. (35) gives the well known result $j_{N}^{(B)}(k=p_{N})=\frac{p_{N}}{\varepsilon_{N}}$ , which reduces to the free current $\frac{p_{N}}{M_{N0}}$ in the nonrelativistic limit.

For the electric current (case $X=E$ in Eq. (34)) we obtain generally

	$\displaystyle j_{b}^{(E)}(k)$	$\displaystyle=v_{b}(k)\,q_{b}+\frac{1}{2}\left(\frac{k}{\varepsilon_{N}}-\frac{\varepsilon_{b}(k)}{\varepsilon_{N}}\,v_{b}(k)\right)$
		$\displaystyle+\frac{p_{N}^{2}}{3\pi^{2}}\,f^{\prime}_{1,bN}(k,p_{N})\,.$		(36)

Here we can insert our model result for $f^{\prime}_{1,bN}$ , given by the last term in Eq. (18) without the factor $\cos\theta$ (see (22)). We can express the result in terms of an effective angular momentum $g$ -factor of the baryon ( $g_{\ell,b}$ ), which we define here — in a naive way — so that it becomes unity for a free proton, i.e., $j_{b}^{(E)}(k\rightarrow 0)\equiv\frac{k}{M_{N0}}\,g_{\ell,b}$ . This gives

\displaystyle g_{\ell,b}=\frac{M_{N0}}{M_{b}}\left(q_{b}-\frac{3x}{1+9x}\left(1+\frac{y_{b}}{2}\right)-\frac{x}{1+x}\,t_{b}\right)\,,

(37)

where $x=2G_{v}\rho_{B}/E_{N}=\frac{1}{9}\left(\frac{\varepsilon_{N}}{E_{N}}-1\right)$ characterizes the strength of the vector interaction. The quantities which depend on the baryon density in (37) are the baryon effective mass $M_{b}$ and $x$ .

III.3 Numerical results

In order to illustrate several physics points of our above discussions, in this subsection we present numerical results for symmetric nuclear matter.

III.3.1 Model parameters

First we explain the choice of our model parameters. The Lagrangian of Eq. (1) contains the coupling constants $G_{\pi}$ and $G_{v}$ , and the current quark masses $m$ and $m_{s}$ , which are related to the constituent quark masses in the vacuum, $M_{0}$ and $M_{s0}$ , by the gap equations (3). The other parameters, which are necessary to define the model, are the infrared (IR) and ultraviolet (UV) cut-offs $\Lambda_{\mathrm{IR}}$ and $\Lambda_{\mathrm{UV}}$ , which are used with the proper-time regularization scheme Schwinger (1951); Hellstern et al. (1997), see App. C. In this scheme, the UV cut-off is necessary to give finite integrals, while the IR cut-off is necessary to avoid unphysical decay thresholds of hadrons into quarks, thereby simulating one important aspect of confinement. These parameters are determined as follows: The IR cut-off should be similar to $\Lambda_{\rm QCD}$ , and we choose $\Lambda_{\rm IR}=0.24\,$ GeV. $\Lambda_{\rm UV}$ , $m$ , and $G_{\pi}$ are determined so as to give a constituent quark mass in vacuum of $M_{0}=0.4\,$ GeV, the pion decay constant $f_{\pi}=0.93\,$ GeV, and the pion mass $m_{\pi}=0.14\,$ GeV, using the standard methods based on the Bethe-Salpeter equation in the pionic $\overline{q}q$ channel Vogl and Weise (1991); Hatsuda and Kunihiro (1994). $m_{s}$ is determined so as to give a constituent $s$ -quark mass in vacuum of $M_{s0}=0.562$ GeV, which reproduces the observed mass of the $\Omega$ baryon $M_{\Omega}=1.67$ GeV by using the quark-diquark bound state equations explained in App. A. The vector coupling $G_{v}$ is determined from the binding energy per-nucleon in symmetric nuclear matter ( $E_{B}/A=-16\,$ MeV) at the saturation density, which becomes $\rho_{0}=0.15\,$ fm^-3. In the present flavor SU(3) NJL model, the vector couplings in the isoscalar and isovector $\overline{q}q$ channels are the same because of constraints from chiral symmetry, and we do not have an independent parameter (like the coupling $G_{\rho}$ in the flavor SU(2) model used in Ref. Tanimoto et al. (2020)) to fit the symmetry energy.⁴⁴4Chiral symmetry would allow different vector couplings in the flavor singlet and octet terms of Eq. (1), but in the mean field approximation used here it is easy to check that there remains only one independent vector coupling in any case. This follows from the identity $\sum_{a=0,3,8}\left(\overline{q}\lambda_{a}\Gamma q\right)^{2}=2\left[\left(\overline{q}_{1}\Gamma q_{1}\right)^{2}+\left(\overline{q}_{2}\Gamma q_{2}\right)^{2}+\left(\overline{q}_{3}\Gamma q_{3}\right)^{2}\right]$ for any Dirac matrix $\Gamma$ . The resulting values of the cut-offs, coupling constants in the $\bar{q}q$ channels, and quark masses are shown in Tab. 1. They are identical to those used in Ref. Tanimoto et al. (2020) except for the $s$ -quark masses which were not needed there. Two additional model parameters are the coupling constants in the scalar and axial vector $qq$ channels, $G_{S}$ and $G_{A}$ of Eq. (49). As explained in App. A, they are fixed to the free nucleon and delta masses ( $M_{N0}=0.94$ GeV, $M_{\Delta 0}=1.23$ GeV). The resulting free masses of octet baryons are then predictions of the model, and are summarized in Tab. 2 together with the observed values.

Table 1: Values for the model parameters which are determined in the vacuum, single hadron, and nuclear matter sectors. The regularization parameters, constituent quark masses in the vacuum (sub-index

0

) and current quark masses are given in units of GeV, and the coupling constants in units of GeV^-2.

$\Lambda_{\rm IR}$	$\Lambda_{\rm UV}$	$G_{\pi}$	$G_{v}$	$M_{0}$	$M_{s0}$	$m$	$m_{s}$
0.240	0.645	19.04	6.03	0.40	0.562	0.016	0.273

Table 2: Masses of octet baryons (in units of GeV) calculated in the vacuum (sub-index

0

) by using the coupling constants

G_{S}=8.76

GeV^-2 and

G_{A}=7.36

GeV^-2 in the scalar and axial vector diquark channels fitted to the vacuum masses of the nucleon and the

\Delta

baryon, in comparison to the observed values.⁶⁶6The fact that our calculated masses agree slightly better with the Gell-Mann Okubo octet mass relation Gell-Mann (1962); Okubo (1962)

M_{N0}+M_{\Xi 0}=\frac{1}{2}\left(M_{\Sigma 0}+3M_{\Lambda 0}\right)

than the experimental values (using either neutral or isospin averaged masses) may be a mere coincidence.

	$M_{N0}$	$M_{\Lambda 0}$	$M_{\Sigma 0}$	$M_{\Xi 0}$
calc.	0.94	1.12	1.17	1.32
obs.	0.94	1.12	1.19	1.32

III.3.2 Energies per nucleon and single baryon energies

In the top panel of Fig. 2 we show the binding energies per nucleon ( ${\cal E}/\rho_{B}-M_{N0}$ ) in symmetric nuclear matter (SNM) in comparison to pure neutron matter (PNM). Although we have only one parameter ( $G_{v}$ ) to fit the binding energy at saturation in SNM, the result for the saturation density agrees with the empirical value. On the other hand, as we do not have any further free parameters, our results for the incompressibility (symmetry energy) in SNM are too large (too small) compared to the empirical values, as will be discussed in more detail in connection to Fig. 3 later. Because of the small symmetry energy, our PNM is slightly bound around densities of $\rho_{B}=0.1$ fm^-3.

In order to show the effects of isospin breaking in PNM on the effective quark and nucleon masses, we list in Tab. 3 the masses in SNM and PNM for four values of the baryon density. Here we can see several points: First, as can be shown from the gap equation (3), for systems with an excess of $d$ -quarks, the magnitude of the mean scalar field $\sigma_{d}$ decreases more rapidly with density than the magnitude of $\sigma_{u}$ . Therefore $M_{u}>M_{d}$ , and one can expect that in an isospin multiplet the baryons with more $u$ -quarks will be heavier. Second, as will be explained in detail later, the isospin splittings for the baryons are generally smaller than for the quarks, because of the scalar isovector polarizabilities of the baryons.

Table 3: Effective masses of quarks and nucleons (in units of GeV) in symmetric nuclear matter (SNM) and pure neutron matter (PNM) for four values of the baryon density (in units of fm^-3).

case	$\rho_{B}$	$M_{u}$	$M_{d}$	$M_{p}$	$M_{n}$
SNM	0	0.4	0.4	0.94	0.94
	0.15	0.325	0.325	0.765	0.756
	0.3	0.284	0.284	0.683	0.683
	0.5	0.257	0.257	0.648	0.648
PNM	0	0.4	0.4	0.94	0.94
	0.15	0.340	0.314	0.768	0.755
	0.3	0.301	0.275	0.695	0.686
	0.5	0.271	0.251	0.656	0.651

In the bottom panel of Fig. 2 we show the Fermi energies (chemical potentials) of the baryons in symmetric nuclear matter. (As the Fermi momenta of hyperons immersed in nuclear matter are zero, the corresponding lines show the energies of hyperons at rest.) The line $\varepsilon_{N}$ shows that nuclear matter is unstable for densities below $0.1$ fm^-3, and at the saturation density it takes the value $940-16=924$ MeV. It is seen that the $\Lambda$ is bound stronger than the nucleon around the saturation density, although its effective mass (not shown here) drops more slowly than $M_{N}$ with increasing density. The reasons are, first, that the curve for $\varepsilon_{\Lambda}$ refers to zero momentum, corresponding to low energy orbitals in finite nuclei, while $\varepsilon_{N}$ refers to the Fermi surface. Second, as shown by Eq. (4), the vector repulsion for the $\Lambda$ in symmetric nuclear matter ( $12\,G_{v}\,\rho_{B}$ ) is smaller than for the nucleon ( $18\,G_{v}\,\rho_{B}$ ), because $\omega^{0}_{s}$ vanishes here.

The curves $\varepsilon_{\Lambda}$ and $\varepsilon_{\Sigma}$ in Fig. 2 show a quite different behavior. Because in this case the vector repulsion is the same, the increase of the difference between the two lines with increasing baryon density reflects the different dependence of their effective masses on $\rho_{B}$ . As a result, around the saturation density the $\Sigma$ is bound by only half of the amount of the $\Lambda$ , i.e., by about 22 MeV less than the $\Lambda$ in our model, which is consistent with the estimate of about 20 MeV presented in Ref. Guichon et al. (2008). The reason for this lies in the different quark substructure: The scalar diquark made of $(u,d)$ quarks, which is the main source of attraction in the $\Lambda$ as well as the nucleon, is absent in the $\Sigma$ as well as in the $\Delta$ baryon. This difference in quark structure, which is well known from the constituent quark model Close (1979), generates the mass difference between the free $\Sigma$ and $\Lambda$ baryons shown in Tab. 2, and increases with increasing baryon density because the mass of the scalar diquark decreases more rapidly than the mass of the axial vector diquark. The strong $(u,d)$ correlations in the scalar channel, as compared to the axial vector channel, play a role similar to the color magnetic spin-spin interaction from gluon exchange. In our model we adjusted this strength to reproduce the $\Delta-N$ mass difference in free space.

The flattening of the curves $\varepsilon_{B}$ with increasing energy, shown in the lower panel of Fig. 2, continues further to the $\Xi$ , because in the present model with 4-fermi interactions the $s$ -quark does not participate in the nuclear interactions in symmetric nuclear matter.

In spite of the increasing $\Sigma-\Lambda$ mass difference due to their different quark substructures, the $\Sigma$ baryon is still bound in our mean field model. It is now believed that the $\Sigma$ is unbound in the nuclear medium Gal et al. (2016), and recent experiments support this view Nanamura et al. (2022). It would be natural as a next step to include the effects of antisymmetrization (exchange terms), both on the level of baryons and the level of quarks. It is, in fact, well known that quark exchange effects appear naturally in the hadronization of the NJL model in the path integral formalism Bentz et al. (2003); Nagata and Hosaka (2004). The effects of the Pauli exclusion principle on the level of quarks to produce the $\Sigma N$ repulsion have been emphasized very much recently Oka (2023); Nanamura et al. (2022). Since the aim of the present work is to explore the effects of the quark substructure of baryons in a mean field approximation for many baryon systems, we will leave this interesting subject for future studies.

III.3.3 Baryon-nucleon Fermi liquid parameters

The top panel of Fig. 3 shows the $\ell=0$ part of the isoscalar baryon-nucleon interaction, given by Eq. (17) without the last term $\propto\bm{p}\cdot\bm{p}^{\prime}$ , and the bottom panel shows the corresponding isovector one, Eq. (18). As in the figure for the baryon energies, the momentum of the nucleons is set to the Fermi momentum ( $p_{N}$ ), and for the hyperons it is set to zero. The behavior of all curves in this figure reflects the change from attraction due to scalar meson exchange at low densities to repulsion from vector meson exchange at higher densities. We find that the third and fourth factors in the first lines of (17) and (18), which reflect the couplings of the scalar mesons to the baryons, decrease substantially in magnitude as the density increases, while the denominators given in the second lines of (17) and (18) become slightly enhanced because of cancellations between the attractive quark loop and repulsive baryon loop contributions. Therefore the attraction from scalar meson exchange decreases much faster with increasing density than for the case of elementary hadrons. In order to illustrate this point more quantitatively, we show in Tab. 4 the various factors which characterize the meson-baryon couplings and meson masses in Eqs. (17) and (18). (The full results for the couplings and meson masses, including the effects of the quark-meson couplings, are given in App. B.)

Table 4: Values of various factors in the effective

\ell=0

interactions

f_{0,BN}

and

f^{\prime}_{0,BN}

of Eq. (17) and (18) which characterize the meson-baryon couplings and meson masses, for four values of the baryon density (in units of fm^-3). The scalar isoscalar polarizability (

\frac{\partial^{2}M_{N}}{\partial M^{2}}

) and the scalar isovector polarizability (

\frac{\partial^{2}M_{p}}{\partial(\Delta M)^{2}}

) of the nucleon are given in units of GeV^-1, and the denominators

D

given in the second lines of (17) and (18) are expressed in units of GeV^-2.

$\rho_{B}$	$\frac{\partial M_{N}}{\partial M}$	$\frac{\partial M_{\Lambda}}{\partial M}$	$\frac{\partial M_{\Sigma}}{\partial M}$	$\frac{\partial M_{\Xi}}{\partial M}$	$\frac{\partial^{2}M_{N}}{\partial M^{2}}$	D
0	2.74	1.83	1.55	0.86	5.90	0.0357
0.15	2.07	1.45	1.15	0.65	12.6	0.0398
0.3	1.49	1.14	0.86	0.49	15.9	0.0532
0.5	1.04	0.89	0.63	0.37	17.2	0.0737
$\rho_{B}$	$\frac{\partial M_{p}}{\partial(\Delta M)}$		$\frac{\partial M_{\Sigma^{+}}}{\partial(\Delta M)}$	$\frac{\partial M_{\Xi^{0}}}{\partial(\Delta M)}$	$\frac{\partial^{2}M_{p}}{\partial(\Delta M)^{2}}$	D
0	0.70		1.55	0.86	9.35	0.0357
0.15	0.49		1.15	0.65	14.7	0.0418
0.3	0.33		0.86	0.49	18.1	0.0572
0.5	0.21		0.63	0.37	19.8	0.0819

The curves in the top panel of Fig. 3 are related to the baryon energies of Fig. 2 by the first of the two general relations given in Eq. (24).⁷⁷7For the case of the nucleon, however, the momentum $k$ is set to the Fermi momentum $p_{N}$ after the differentiation in (24). It is thus natural that the average values of $f_{0,BN}$ become smaller in the sequence $NN\rightarrow\Lambda N\rightarrow\Sigma N\rightarrow\Xi N$ . In particular, as explained above, $M_{\Lambda}$ decreases faster with density than $M_{\Sigma}$ , and therefore the first term in (28) shows that the $\Lambda N$ attraction at low densities is stronger than the $\Sigma N$ attraction. Around the saturation density, the $\Lambda N$ and the isoscalar $\Xi N$ interactions are similar and very small, while the isoscalar $NN$ and $\Sigma N$ interactions are both repulsive.⁸⁸8We remind again that $f_{0,BN}$ refers to the spin averaged interaction characterized by $\ell=0$ and $T=0$ in the particle-hole channel. For the $NN$ case, we can use Eqs. (25) and (28) to split the incompressibility as $K=\frac{3p_{N}^{2}}{E_{N}}+9\rho_{B}\frac{M_{N}}{E_{N}}\frac{\partial M_{N}}{\partial\rho_{B}}+162\,G_{v}\,\rho_{B}=(0.253-1.014+1.124)$ GeV = $0.363$ GeV, where the first term refers to noninteracting quasiparticles with $E_{N}=0.8$ GeV, the second term corresponds to $\sigma$ meson exchange and the third term to $\omega$ meson exchange. In order to reproduce the empirical value $K\simeq 0.25$ GeV, we would need $f_{0,NN}\simeq 0$ at saturation density, instead of the positive value indicated in Fig. 3.

The 3 curves in the lower panel of Fig. 3 similarly result from the attraction due to $\delta$ meson exchange at low densities and the repulsion from $\rho$ meson exchange at higher densities. The fact that the isovector $\Sigma N$ repulsion is stronger than the others is simply because of the isospin factor $t_{b}$ in Eq. (29), which indicates that the energy of $\Sigma^{\pm}$ is most sensitive to changes of the isovector nucleon density. The overall size of the isovector interactions is small compared to the isoscalar ones. For the $NN$ case, we can use Eqs. (26) and (29) to split the symmetry energy as $a_{s}=\frac{p_{N}^{2}}{6\,E_{N}}+\frac{\rho_{B}}{2}\frac{M_{N}}{E_{N}}\frac{\partial M_{p}}{\partial\rho_{(3)}}+G_{v}\,\rho_{B}=(14-3+7)$ MeV = $18$ MeV, where the first term refers to noninteracting quasiparticles, the second term corresponds to $\delta$ meson exchange and the third term to $\rho$ meson exchange. It is known from the case of elementary nucleons Ulrych and Muther (1997) that the mechanism of $\delta$ meson exchange gives a negative contribution to the symmetry energy, and in our model this effect is small. Our value of $a_{s}$ is considerably smaller than the empirical value $a_{s}\simeq 32$ MeV, which reflects the fact that our 3-flavor Lagrangian (1) does not allow for an independent vector coupling in the isovector channel because of the assumed flavor and chiral symmetry, in contrast to the 2-flavor case Tanimoto et al. (2020).

Finally in this subsection, we add two more comments. The first concerns the isospin splittings which can be expected for isospin asymmetric matter. Because our $f^{\prime}_{0,BN}$ is negative at small densities, the first term in Eq. (29) is negative for $b=p,\Sigma^{+},\Xi^{0}$ . For systems with neutron excess ( $\rho_{(3)}<0$ ) we can then expect that the in-medium isospin splittings will be ordered such that the particles with more $u$ -quarks become heavier, which is consistent with our finding that $M_{p}>M_{n}$ and $M_{u}>M_{d}$ in neutron rich matter, see Tab. 3. The reason why the mass splittings for baryons are smaller than for quarks is now clear from Tab. 4, which shows that the isovector couplings $\partial M_{b}/\partial(\Delta M)$ strongly decrease with increasing baryon density. Expressed in a different way, the scalar isovector polarizability of the nucleon ( $\frac{\partial^{2}M_{p}}{\partial(\Delta M)^{2}}$ ) strongly increases with the density.

Second, it is well known that any two-body interaction with non-explicit density dependence, for example through masses and couplings, contains effects from an effective three-body interaction. Taking the $\ell=0$ part of Eq. (17) as an example, in the case of point nucleons the only density dependence of this kind resides in the factor $\frac{M_{N}}{E_{N}}$ and in the function $\phi_{N}$ in the denominator. The decrease of our couplings and the increase of meson masses due to the scalar isoscalar polarizability of the nucleons reflect the presence of additional repulsive three-body interactions.⁹⁹9The variation of $f_{\,0,bN}$ with density can be expressed as an effective three-body interaction: $\displaystyle\frac{\delta f_{0,bN}}{\delta\rho_{B}}=\frac{1}{4}\sum_{\tau=p,n}\left(h_{0,bp\tau}+h_{0,bn\tau}\right)\,,$ where the three-particle forward scattering amplitudes ( $h_{0}$ ) are defined as averages over the angles between the momenta of the three interacting particles Bentz and Cloët (2022). The rapid decrease of the $bN$ attraction with increasing density, expressed by Fig. 3, shows that our effective three-particle interaction is strongly repulsive, but — as Fig. 2 (lower panel) has shown — not sufficient to generate an overall repulsion between the $\Sigma$ baryon and the nucleon.

III.3.4 In-medium orbital $g$ -factors of baryons

Here we wish to illustrate the renormalization of the orbital angular momentum $g$ -factors, given by Eq. (37), for a few cases. In Sec. III.2.3 we used the concept of the backflow, which is central to the Fermi liquid theory, but the same results can be obtained in relativistic meson-nucleon theories by using the response of the core (filled Fermi sea of nucleons) to the addition of one nucleon McNeil et al. (1986); Ichii et al. (1988) or one hyperon Cohen and Furnstahl (1987); Cohen (1993). In such a description, the backflow arises from RPA-type vertex corrections due to virtual $N\overline{N}$ excitations of the core,¹⁰¹⁰10The “antinucleons” which show up in those vertex corrections, or in the Z-graph contributions to the scalar meson propagators mentioned in Sec. III.1, are highly virtual objects, mathematically necessary to form a complete set of spinors, and have little to do with real observable antinucleons. and the importance of these contributions to give reasonable magnetic moments in relativistic theories is well known Arima et al. (1987); Furnstahl and Serot (1987). As examples for baryons with positive charge, we illustrate the relation (37) for the proton and the $\Sigma^{+}$ , and as examples for neutral baryons we show the cases of the neutron and the $\Lambda$ in Fig. 4. For the isoscalar combination $g_{\ell,p}+g_{\ell,n}$ , the backflow reduces the enhancement ( $\frac{M_{N0}}{M_{N}}\simeq 1.24$ near the saturation density) by a factor of $\frac{E_{N}}{\varepsilon_{N}}\simeq 0.87$ , while for the isovector combination $g_{\ell,p}-g_{\ell,n}$ there is almost no reduction, because the last term in (37) is very small. As a result, the isovector combination remains enhanced, i.e., $g_{\ell}$ of the proton (neutron) is larger (smaller) than its free value. For the $\Sigma^{+}$ , the enhancement due to its reduced mass is only about half of the case of the proton, and the reduction from the backflow gives results which change only mildly with density. For the $\Lambda$ , the backflow corrections are similar in magnitude to the case of the neutron, but its effective mass, and therefore also $g_{\ell}$ , decreases more slowly with density. For more expensive discussions on backflow effects for the magnetic moments of hypernuclei, we refer to Ref. Cohen (1993).

III.3.5 Comments on sizes of quark cores of in-medium nucleons

Finally in this section, we wish address the question whether the size of in-medium nucleons invalidates the basic physical picture of the mean field approximation. The relevance of this question is underlined by the fact that the NJL model is known to predict a moderate swelling of nucleons in the medium at normal densities, a feature which is important for the EMC effect Cloet et al. (2006) or the Coulomb sum rule Cloët et al. (2016). If the nucleons swell considerably at higher densities, the Pauli principle would become inapplicable at the nucleon level.

Rather than the physical size of nucleons including their meson clouds, the quantity which seems more relevant for role of the Pauli principle is the size of the quark cores of the nucleons in the medium. Here we consider the rms radius of the baryon density distribution of the quark cores, denoted as $r_{N}(\rho_{B})$ , which is an isoscalar quantity and therefore the same for protons and neutrons. The definitions and further details are given in App. D, and the results are shown in Fig. 5. Our free nucleon (zero density) value is $r_{N}(0)=0.475$ fm, which increases by $7\%$ at saturation density ( $0.15\,$ fm^-3), and by $13\%$ at $\rho_{B}=0.5\,$ fm^-3. Even at very large densities ( $\rho_{B}\simeq 1.0$ fm^-3) the baryon radius of the quark core increases only by $16\%$ of its free value. This behavior reflects our phenomenological implementation of confinement effects via the infrared cut-off ( $\Lambda_{\rm IR}$ ). It is interesting to note that our values of $r_{N}$ are similar to the radii which have been assumed in the excluded volume framework in QMC model calculations Panda et al. (2002); Leong et al. (2023a), although we do not go into further details here.

By using the values of $r_{N}$ shown in Fig. 5, we can estimate the volume fractions occupied by the quark cores in SNM (see App. D). We obtain $9\%$ at saturation density, and $36\%$ at the highest density shown in the previous figures of this section ( $\rho_{B}=0.5$ fm^-3). Although these numbers may give us some confidence in the overall physical picture of the mean field approximation, they leave room for corrections and improvements of the model. We also remind that the Pauli principle at the quark level has been predicted to play an important role in producing the $\Sigma N$ repulsion even at normal densities, as mentioned at the end of Sec. III.3.2. Further investigations on these points are necessary.

IV NEUTRON STAR MATTER

In this section we wish to discuss our results for neutron star matter and the resulting star masses, based on the expression (6) for the energy density and the equilibrium and charge neutrality conditions (11). Our parameters are the same as used in symmetric nuclear matter, see Sec. III.3.1. As mention at the end of Sec. II.2, we will not analyze the effective baryon-baryon interactions in neutron star matter as exhaustively as we have done it for nuclear matter, in order to keep the length of the paper within reasonable limits.

IV.1 Single particle properties in-medium

First we show our results for the quark effective masses in Fig. 6 as functions of the baryon density. Because of the isospin asymmetry (excess of $d$ -quarks) in neutron star matter, the $u$ -quark becomes heavier than the $d$ -quark by 25 MeV at baryon densities around $0.3$ fm^-3. As discussed already in Sec. III.3.2, this is expected from $|\sigma_{u}|>|\sigma_{d}|$ in neutron rich matter, or equivalently from the effective $\delta$ -meson exchange mechanism Ulrych and Muther (1997) in hadronic theories. The $s$ -quark mass, on the other hand, starts to decrease as soon as hyperons appear in the system, i.e., as soon as the condition $\partial{\cal E}/\partial\sigma_{s}=0$ receives contributions from hyperons in the baryon loop term ${\cal E}_{B}$ of Eq. (6). In this case, the $\bar{s}s$ exchange between hyperons can proceed without violating the OZI rule Okubo (1963); Zweig (1964); Iizuka (1966), and, as anticipated in Sec. I, this gives rise to an appreciable attraction in neutron star matter. We will explain later how this decrease of $M_{s}$ influences the masses of neutron stars.

The upper panel of Fig. 7 shows our results for the chemical potentials, and the lower panel shows the density fractions of the particles, as functions of the baryon density. The upper panel is the analogue of the lower panel of Fig. 2, discussed in the previous section for symmetric nuclear matter. The three solid lines in the upper panel of Fig. 7 (from bottom to top) show $\mu_{p}=\mu_{n}-\mu_{e}$ , $\mu_{n}$ , and $\mu_{n}+\mu_{e}$ . Because of the conditions (11), the density where the chemical potential of a hyperon with electric charge $q_{b}$ touches the solid line $\mu_{n}-q_{b}\mu_{e}$ from above is the threshold density for this hyperon. Below the threshold densities, the chemical potentials are simply the energies of hyperons at rest (zero Fermi momentum). Compared to the symmetric nuclear matter case of Fig. 2, the lines show a considerable isospin splitting, which mainly comes from the vector potential term in Eq. (4). For example, the vector potential for the $\Sigma^{-}$ is $4G_{v}\left(2\rho_{d}+\rho_{s}\right)$ , which is larger than the vector potential for the $\Sigma^{+}$ , which is $4G_{v}\left(2\rho_{u}+\rho_{s}\right)$ . The mass splittings are in the opposite order, e.g., the $\Sigma^{+}$ is heavier than the $\Sigma^{-}$ , as can be expected also from Fig. 6. The mass splittings in baryon isospin multiplets are, however, small compared to the splittings from the vector potential. For example, at baryon densities around $0.3$ fm^-3 the mass splitting between the $\Sigma^{+}$ and the $\Sigma^{-}$ is only about $20$ MeV, and the proton-neutron mass difference is only about $10$ MeV, both being smaller than the naive expectation from the quark mass difference shown in Fig. 6 for the reasons explained in the previous section. In the low density region, where $\rho_{s}=0$ and $\rho_{u}+\rho_{d}=3\rho_{B}$ , the vector potentials for the $\Sigma^{0}$ and the $\Lambda$ are the same ( $4G_{v}(\rho_{u}+\rho_{d})$ ), and we see again the different behaviors of their energies with increasing baryon density, which is caused by their different quark substructures as discussed in Sec. III.3.2.

As we can see from Fig. 7, the threshold density for the $\Sigma^{-}$ is $\rho_{B}=0.35$ fm^-3 in our calculation. Although it has been conjectured for long on energetic reasons that the $\Sigma^{-}$ will appear as the first hyperon in neutron star matter Heiselberg and Hjorth-Jensen (2000), this point is controversial nowadays Stone et al. (2021); Motta and Thomas (2022); Leong et al. (2023b), mainly because the $\Sigma N$ interaction is believed to be repulsive (see the related discussions at the end of Sec. III.3.2). However, we wish to note that also in our present mean field model the onset of the $\Sigma^{-}$ depends on several details: First, we are underestimating the free $\Sigma$ mass by about 20 MeV (see Tab. 2); second, the in-medium mass of the $\Sigma^{-}$ is shifted down by a similar amount relative to the $\Sigma^{+}$ as explained above; and third, our electron chemical potential is rather large in this density region. Therefore, apart from the more fundamental problem on the $\Sigma N$ repulsion, the question whether $\mu_{\Sigma^{-}}$ touches $\mu_{n}+\mu_{e}$ or not, and if it does at which baryon density, depends on several details of the model. (We will return to this point in a different context in Sec. VI.)

IV.2 Equation of state and neutron star masses

The upper panel of Fig. 8 shows our results for the pressure in neutron star matter as function of the baryon density, and the lower panel shows the neutron star masses, as obtained from the solution of the Tolman-Oppenheimer-Volkoff (TOV) equations Tolman (1939); Oppenheimer and Volkoff (1939), with the constraints of Eq. (11) imposed, as functions of the central baron density. We show the cases of nucleons and leptons only, the case including hyperons but fixing artificially the effective $s$ -quark mass to its free value ( $M_{s0}$ ), and the full result with the $s$ -quark mass determined by minimization of the energy density. The results for nucleons and leptons only are very similar to the results obtained in Ref. Tanimoto et al. (2020) for the flavor $SU(2)$ case, although there it was possible to reproduce the symmetry energy without explicit breaking of chiral symmetry of the interaction Lagrangian. (See App. A for a more detailed comparison.)

It is well known that the presence of hyperons can lead to a sizable reduction of the pressure in neutron star matter and a decrease of the maximum mass of neutron stars Glendenning (1997); Bombaci (2017), and Fig. 8 shows that the same situation is encountered in a relativistic mean field calculation which takes into account the internal quark-diquark structure of the octet baryons. Our results suggest that most of the reduction of the pressure arises simply because nucleons and leptons with high Fermi momenta can be converted to hyperons with low Fermi momenta by weak processes. The reduction of the $s$ -quark mass in the medium is not so important for the overall size of the pressure and the maximum star mass, but it works towards destabilization of the star as the central baryon density increases. The values of the maximum central baryon densities which gives stable stars, the maximum star masses, and the radii of the stars with maximum mass for the three cases shown in Fig. 8 are as follows:

\displaystyle(\rho_{B}^{\rm max}(r=0),\,M_{\rm star}^{\rm max},\,R)=(0.9\,{\rm fm}^{-3},\,\,2.17\,M_{\odot},\,\,11.5\,{\rm km})

for the case of no hyperons;

\displaystyle(\rho_{B}^{\rm max}(r=0),\,M_{\rm star}^{\rm max},\,R)=(0.85\,{\rm fm}^{-3},\,\,1.83\,M_{\odot},\,\,11.8\,{\rm km})

for the case with hyperons but $M_{s}$ fixed to $M_{s0}$ ; and

\displaystyle(\rho_{B}^{\rm max}(r=0),\,M_{\rm star}^{\rm max},\,R)=(0.72\,{\rm fm}^{-3},\,\,1.73\,M_{\odot},\,\,12.3\,{\rm km})

for the case with hyperons and $M_{s}$ determined from minimization of the energy density.

V ROLE OF 6-FERMI AND 8-FERMI INTERACTIONS

Because the maximum mass of neutron stars is sensitive to the high density behavior of the equation of state, it is natural to investigate the role of higher order Fermi interactions, i.e., the 6-fermi ’t Hooft (1976) and 8-fermi Osipov et al. (2007) interactions. While there is little doubt about the importance of the 6-fermi (flavor determinant) interaction to break the $U_{A}(1)$ symmetry of the 4-fermi Lagrangian of Eq. (1) and to split the masses of the otherwise degenerate pseudoscalar mesons $\pi$ and $\eta$ , the situation is not so clear for the 8-fermi interactions, because many possible flavor structures are allowed by chiral symmetry. In this work we will limit ourselves to three types of chiral invariant 8-fermi interactions with the simplest structure, namely the square of the scalar-pseudoscalar term in Eq. (1), the product of this term with the vector-axial vector term, and the square of the vector-axial vector term. We wish to investigate whether those higher order Fermi interactions in the $\overline{q}q$ channels, with coupling constants restricted by the basic properties of symmetric nuclear matter around the saturation point, can lead to appreciable changes in high density neutron star matter or not. We will not include higher order interactions in the $qq$ channels used to construct the baryons as quark-diquark bound states, i.e., the Lagrangian of Eq. (49) is left unchanged.

V.1 Basic formulas and new parameters

To the basic NJL Lagrangian of Eq. (1), we add the 6-fermi (flavor determinant) interaction Buballa (2005)

\displaystyle{\cal L}_{6}

\displaystyle=G_{6}\,{\rm det}\left[\bar{q}_{\alpha}\left(1-\gamma_{5}\right)q_{\beta}+\bar{q}_{\alpha}\left(1+\gamma_{5}\right)q_{\beta}\right]\,,

(38)

and the following 8-fermi interactions:

	$\displaystyle{\cal L}_{8}$	$\displaystyle=G^{(ss)}_{8}\,\left({\cal L}_{s}\,{\cal L}_{s}\right)-G^{(sv)}_{8}\,\left({\cal L}_{s}\,{\cal L}_{v}\right)-G^{(vv)}_{8}\,\left({\cal L}_{v}\,{\cal L}_{v}\right)$
		$\displaystyle\equiv{\cal L}_{8}^{(ss)}+{\cal L}_{8}^{(sv)}+{\cal L}_{8}^{(vv)}\,.$		(39)

Here ${\cal L}_{s}$ means the 4-fermi interaction in the scalar-pseudoscalar channel of Eq. (1) without the factor $G_{\pi}$ , and ${\cal L}_{v}$ means the one in the vector-axial vector channel of Eq. (1) without the factor $(-G_{v})$ . In this simplest possible form, each factor ${\cal L}_{s}$ or ${\cal L}_{v}$ is closed under the summations over Dirac, flavor, and color indices. Altogether 4 new coupling constants are involved in (38) and (39).

The mean field approximation is implemented in the same way as for the 4-fermi interactions in Sec. II.1. The gap equation (3) is now replaced by the more complicated form

\displaystyle M_{\alpha}=m_{\alpha}-\sigma_{\alpha}\left(1+\frac{G_{8}^{(ss)}}{4G_{\pi}^{3}}\,\sigma_{\beta}^{2}-\frac{G_{8}^{(sv)}}{8G_{\pi}G_{v}^{2}}\,\omega_{\beta}^{2}\right)+\frac{G_{6}}{8G_{\pi}^{2}}\,\sigma_{\beta}\,\sigma_{\gamma}\,,

(40)

where in the 6-fermi term $(\alpha,\beta,\gamma)$ is any set of three different quark flavors, and in the other terms a sum over the quark flavors $\beta$ is implied. The baryon energies (4) are replaced by

\displaystyle\varepsilon_{b}(k)=E_{b}(k_{b})+n_{\alpha/b}V_{\alpha}^{0}\,,

(41)

where $\bm{k}_{b}=\bm{k}-n_{\alpha/b}\bm{V}_{\alpha}$ , with the vector fields $V_{\alpha}^{\mu}$ defined by

\displaystyle V_{\alpha}^{\mu}=\omega_{\alpha}^{\mu}\left(1+\frac{G_{8}^{(sv)}}{8G_{\pi}^{2}G_{v}}\,\sigma_{\beta}^{2}+\frac{G_{8}^{(vv)}}{4G_{v}^{3}}\,\omega_{\beta}^{2}\right)\,.

(42)

The new contributions from the 6-fermi and 8-fermi interactions to the energy density are

$\displaystyle{\cal E}_{6}$	$\displaystyle=-\frac{G_{6}}{16G_{\pi}^{3}}\,\left(\sigma_{u}\,\sigma_{d}\,\sigma_{s}-\sigma_{u0}\,\sigma_{d0}\,\sigma_{s0}\right)\,,$	(43)
$\displaystyle{\cal E}_{8}$	$\displaystyle=\frac{3G_{8}^{(ss)}}{64G_{\pi}^{4}}\,\left(\sigma_{\alpha}^{2}\,\sigma_{\beta}^{2}-\sigma_{\alpha 0}^{2}\,\sigma_{\beta 0}^{2}\right)-\frac{3G_{8}^{(sv)}}{64G_{\pi}^{2}G_{v}^{2}}\,\sigma_{\alpha}^{2}\,\omega_{\beta}^{2}$
	$\displaystyle-\frac{3G_{8}^{(vv)}}{64G_{v}^{4}}\,\omega_{\alpha}^{2}\,\omega_{\beta}^{2}\,,$	(44)

which are added to Eq. (6), after replacing $\omega_{\alpha}^{\mu}$ in (6) by the expression given in Eq. (42). It is easy to check that the basic conditions (8), which determine the three scalar and three vector mean fields $\sigma_{\alpha}$ and $\omega_{\alpha}^{\mu}$ , lead to the same expressions (9) and (10) as before, because those expressions simply reflect the definitions given by Eq. (2). If we eliminate the vector fields by using (10), it becomes clear that $G_{8}^{(ss)}$ and $G_{8}^{(vv)}$ must be positive in order that the energy density is bounded from below, while the sign of $G_{8}^{(sv)}$ is not determined generally. For the case of neutron star matter, the conditions of chemical equilibrium and charge neutrality are given by Eq. (11) with the modified baryon chemical potentials $\mu_{b}=\varepsilon_{b}(k=p_{b})$ .

We also note that the 6-fermi and 8-fermi interactions lead to a renormalization of the residual 4-fermi interactions. The only physical quantities, for which we use the residual 4-fermi interactions in the $\overline{q}q$ channel to fix model parameters in this work, are the mass of the pion, the pion decay constant, and the $\eta$ - $\eta^{\prime}$ mass difference, where the pseudoscalar mesons $\eta$ and $\eta^{\prime}$ arise from mixing Hatsuda and Kunihiro (1991); Rehberg et al. (1996) between the $\eta_{0}$ and $\eta_{8}$ . The effective 4-fermi coupling constants in the vacuum, relevant for those quantities, are given by (see, for example, Refs. Kato et al. (1993); Rehberg et al. (1996) for the 6-fermi case)

$\displaystyle\tilde{G}_{\pi}$	$\displaystyle=\left(G_{\pi}+\frac{G_{8}^{(ss)}}{4G_{\pi}^{2}}\,\sigma_{\alpha 0}^{2}\right)-\frac{G_{6}}{8G_{\pi}}\,\sigma_{s0}\equiv 19.04\,{\rm GeV}^{-2}\,,$	(45)
$\displaystyle\tilde{G}_{00}$	$\displaystyle=\left(G_{\pi}+\frac{G_{8}^{(ss)}}{4G_{\pi}^{2}}\,\sigma_{\alpha 0}^{2}\right)+\frac{G_{6}}{12G_{\pi}}\left(2\sigma_{0}+\sigma_{s0}\right)\,,$	(46)
$\displaystyle\tilde{G}_{08}$	$\displaystyle=-\frac{\sqrt{2}\,G_{6}}{12G_{\pi}}\left(\sigma_{0}-\sigma_{s0}\right)\,,$	(47)
$\displaystyle\tilde{G}_{88}$	$\displaystyle=\left(G_{\pi}+\frac{G_{8}^{(ss)}}{4G_{\pi}^{2}}\,\sigma_{\alpha 0}^{2}\right)-\frac{G_{6}}{24G_{\pi}}\left(4\sigma_{0}-\sigma_{s0}\right)\,.$	(48)

We require that $\tilde{G}_{\pi}$ has the same value as $G_{\pi}$ in the 4-fermi calculation in order to reproduce the observed pion mass (see Tab. 1), and that $G_{6}$ reproduces the observed mass difference $m_{\eta^{\prime}}-m_{\eta}=0.41$ GeV. One can use Eq. (45) to express the quantity $\left(G_{\pi}+\frac{G_{8}^{(ss)}}{4G_{\pi}^{2}}\,\sigma_{\alpha 0}^{2}\right)$ in the form $\tilde{G}_{\pi}+\frac{G_{6}}{8G_{\pi}}\,\sigma_{s0}$ . By inserting this into (46) and (48), we see that the three coupling constants (46), (47) and (48), which are used to calculate the $\eta-\eta^{\prime}$ mass difference, can be expressed in terms of $\tilde{G}_{\pi}$ , $G_{6}$ , and the quark condensates in the vacuum, $\langle\overline{u}u\rangle_{0}=\langle\overline{d}d\rangle_{0}=\sigma_{0}/(4G_{\pi})$ and $\langle\overline{s}s\rangle_{0}=\sigma_{s0}/(4G_{\pi})$ , which are fixed by the constituent quark masses in the vacuum and the cut-offs given in Tab. 1. Therefore $G_{6}$ can be adjusted to the $\eta-\eta^{\prime}$ mass difference in the standard way Kato et al. (1993), without recourse to the value assumed for $G_{8}^{(ss)}$ . It is also easy to see that the gap equation (40) for the $u,d$ quarks in the vacuum remains numerically the same as in the pure 4-fermi case, because it can be expressed as $M_{0}=m-4\tilde{G}_{\pi}\langle\overline{u}u\rangle_{0}$ . Therefore the value of $m$ , given in Tab. 1, is unchanged.¹¹¹¹11The value of the current $s$ -quark mass depends slightly on the values assumed for $G_{6}$ and $G_{8}^{(ss)}$ . Also the original 4-fermi coupling constant $G_{\pi}$ changes according to Eq. (45), although this has no effect on any physical quantity. By the standard calculations, we find that $G_{6}=1260$ GeV^-5 reproduces the observed $\eta-\eta^{\prime}$ mass difference.

Next we comment on the role of the 8-fermi coupling constants. As one can expect from the gap equation (40), $G_{8}^{(ss)}$ works into the same direction as the original 4-fermi coupling $G_{\pi}$ , i.e., it gives attraction, while a positive coupling $G_{8}^{(sv)}$ gives repulsion. The coupling $G_{8}^{(vv)}$ , on the other hand, is not related to the gap equation, but after eliminating the vector fields according to (10) it is easily seen to give a repulsive contribution of $4G_{8}^{(vv)}\rho_{\alpha}^{2}\rho_{\beta}^{2}=81G_{8}^{(vv)}\rho_{B}^{4}$ to the energy density, and $12\,G_{8}^{(vv)}\rho_{\alpha}^{2}\rho_{\beta}^{2}=243\,G_{8}^{(vv)}\rho_{B}^{4}$ to the pressure in symmetric nuclear matter. Although the 8-fermi coupling constants can be treated as free parameters, their choice is strongly limited by the requirements that the saturation point of isospin symmetric nuclear matter is unchanged, and the discrepancies of the calculated incompressibility and the symmetry energy to the empirical values do not increase. In the present calculation, we achieved this by making use of the balance between the attractive $(ss)$ -type interaction and the repulsive $(vv)$ -type interaction. Concerning the $(sv)$ -type interaction, which can work as an attraction ( $G_{8}^{(sv)}<0$ ) or a repulsion ( $G_{8}^{(sv)}>0$ ), we found that the case of attraction leads to conflicts with the nuclear matter equation of state, and the case of repulsion gives a much smaller effect in neutron star matter than the repulsive $(vv)$ -type interaction. We therefore consider only the case $G_{8}^{(sv)}=0$ in the calculations described below. We also note that changes in the original 4-fermi vector coupling constant $G_{v}$ , under the constraints imposed by symmetric nuclear matter, do not lead to any noteworthy improvements of the equation of state of neutron star matter, so we keep the same value as given in Tab. 1.

In Tab. 5, we list as case 1 the pure 4-fermi case, where the 6-fermi and 8-fermi coupling constants are zero, and in case 2 the 6-fermi interaction with the value of $G_{6}$ as determined above is added. Case 3 gives one possible choice for the 8-fermi coupling constants, where the balance between the attractive $(ss)$ -type interaction and the repulsive $(vv)$ -type interaction is used to keep the nuclear matter properties around the saturation point unchanged, while the $(sv)$ -type interaction is assumed to vanish.

Table 5: Values of the 6-fermi coupling constant

G_{6}

in units of GeV^-5, and the 8-fermi coupling constants

G_{8}^{(ss)}

and

G_{8}^{(vv)}

in units of GeV^-8, for the three cases discussed in this section. The coupling

G_{8}^{(sv)}

is set to zero in all three cases.

case	$G_{6}$	$G_{8}^{(ss)}$	$G_{8}^{(vv)}$
1	0	0	0
2	1260	0	0
3	1260	2330	1220

V.2 Numerical results including 6-fermi and 8-fermi interactions

The top panel of Fig. 9 shows the results for the pressure in neutron star matter for the three cases listed in Tab. 5 as functions of the baryon density, and the bottom panel shows the resulting neutron star masses as functions of the central baryon density. Although the 6-fermi interaction (case 2) leads only to a slight decrease of the pressure in the region of $\rho_{B}=0.35\sim 0.8$ fm^-3, the resulting decrease in neutron star masses is quite significant.¹²¹²12As in other works, for example Ref. Weissenborn et al. (2012), we find that small changes of the pressure in this region of baryon densities can lead to appreciable changes in the star masses. On the other hand, the $(vv)$ -type 8-fermi interaction (last term in Eq. (39)), with a very moderate coupling constant and counterbalanced by the $(ss)$ -type interaction so as not to change the saturation properties of symmetric nuclear matter, gives a strongly increasing pressure for $\rho_{B}>0.7$ fm^-3 and stabilizes the neutron stars against collapse for central densities larger than $0.7$ fm^-3.

Taken together, the 6-fermi and 8-fermi interactions do not lead to noticeable changes of the maximum star masses, but rather work towards stabilization of the stars for high central densities. The resulting star masses for case 3 in the range of central densities between $0.7$ and $1.5$ fm^-3 are all around $1.7\,M_{\odot}$ , with radii decreasing from $11.5$ km to $9.5$ km. We finally give the values of the maximum central baryon densities which gives stable stars, the maximum star masses, and the radii of the stars with maximum mass for the three cases shown in Fig. 9:

\displaystyle(\rho_{B}^{\rm max}(r=0),\,M_{\rm star}^{\rm max},\,R)=(0.72\,{\rm fm}^{-3},\,\,1.73\,M_{\odot},\,\,12.3\,{\rm km})

for the case 1;

\displaystyle(\rho_{B}^{\rm max}(r=0),\,M_{\rm star}^{\rm max},\,R)=(0.8\,{\rm fm}^{-3},\,\,1.62\,M_{\odot},\,\,11.9\,{\rm km})

for the case 2; and

\displaystyle(\rho_{B}^{\rm max}(r=0),\,M_{\rm star}^{\rm max},\,R)=(1.4\,{\rm fm}^{-3},\,\,1.72\,M_{\odot},\,\,9.8\,{\rm km})

for the case 3.

VI SUMMARY

In this paper we used the 3-flavor NJL model as an effective quark theory of QCD to describe the octet baryons as quark-diquark bound states, and the equations of state of nuclear and neutron star matter in the relativistic mean field approximation based on quark degrees of freedom. One of our basic concepts was to preserve the flavor and chiral symmetries of the interaction Lagrangian, i.e., to allow explicit symmetry breakings only by the current quark masses and not by ad-hoc changes of model parameters. In Sec. I we stated the four main purposes of our work, so let us now summarize our results in this order.

First, the internal quark structure of baryons leads to density dependent meson-baryon coupling constants and meson masses which strongly reduce the attractive parts of the interactions in nuclear matter. The main reason for this effect is the nonlinear behavior of the hadron masses as functions of the constituent quark masses. In particular, we found that the attraction experienced by the $\Sigma$ baryon immersed in nuclear matter is reduced more strongly than that for the $\Lambda$ baryon, and we could verify that the mass difference between the $\Sigma$ and the $\Lambda$ baryons immersed in the nuclear medium increases with increasing density. However, we found that this effect, which is based on the different quark-diquark structures of those two baryons, is not sufficient to make the $\Sigma$ unbound in the region of normal nuclear matter density.

Second, we used concepts of the relativistic Fermi liquid theory to derive the effective meson exchange interaction between octet baryons in the nuclear medium, and the analogue of the Landau relation between the energies of the baryons and the interactions between them. We also used the same concepts to discuss the renormalization of currents carried by baryons, as well as the effects of nucleon density variations on the energies of hyperons immersed in the nuclear medium. To the best of our knowledge, some of these relations cannot be found in the literature, and we hope that our results will be useful for further investigations.

Third, we designed our mean field approximation so that it reflects the basic symmetries of the model and their dynamical breakings, regardless of possible disagreements with observations. To appreciate this point, let us suppose for the moment that we had explicitly broken the flavor and chiral symmetries, as specified below Eq. (1), by choosing a different coupling constant (say $G_{\rho}$ ) for the isovector term in the second line of Eq. (1): $G_{\rho}\,\left[\left(\overline{q}\lambda_{i}\,\gamma^{\mu}\,q\right)^{2}+\left(\overline{q}\lambda_{i}\,\gamma^{\mu}\gamma_{5}\,q\right)^{2}\right]$ , where $i=1,2,3$ . By choosing $G_{\rho}\simeq 3G_{v}\simeq 18$ GeV^-2, we could reproduce the empirical symmetry energy $a_{s}=32$ MeV (see Sec. III.3.3), the shallow bound state of pure neutron matter in Fig. 2 would disappear, neutron stars made of nucleons and leptons would become heavier, and the onset of the $\Sigma^{-}$ baryon would move to higher densities or disappear, because its energy gets a positive shift from the vector isovector potential, twice as large as for the neutron (see upper panel of Fig. 7, and Eq.(27)). This would delay the onset of the decrease of $M_{s}$ in neutron star matter (see Fig. 6) and thereby hinder the succession of further hyperons (see lower panel of Fig. 7), leading again to larger star masses. While this ad-hoc modification may still have some phenomenological justification, one may think of more drastic changes, like for example enhancing the coupling constant in the vector potential acting on the $s$ -quark in Eq. (2), or introducing a phenomenological repulsive function into the energy density which grows asymptotically for large densities. In these or other ways one could “improve” the results, but only little can be learned from it.

Fourth, we found that the so called hyperon puzzle persists in the NJL model for composite octet baryons in the mean field approximation, and 6-fermi and 8-fermi interactions - with coupling constants chosen so as not to spoil the saturation properties of normal nuclear matter - do not solve the problem. On the positive side, we have shown that a special kind of 8-fermi interaction, characterized by a product of four quark current operators, is able to support stable stars up to $1.7$ solar masses over a large region of central densities. In view of the extremely large baryons densities involved in the investigation of neutron stars, we believe that any solution to the hyperon puzzle must involve quark degrees of freedom, not only quarks in individual hadrons but also quarks which belong to two or more hadrons, or to the whole system. An investigation along these lines would naturally lead to an examination of various patterns of phase transitions to 3-flavor quark matter, including pairing and condensation phenomena.

Acknowledgements.

K. N. wishes to thank the staff and students of the Department of Physics at Tokai University for their discussions and advice. W.B. acknowledges very helpful advice from Prof. H. Tamura and Prof. F. Weber. The work of I.C. was supported by the U.S. Department of Energy, Office of Science, Office of Nuclear Physics, contract no. DE-AC02-06CH11357. This work was supported partially through USJHPE (U.S. - Japan Hadronic Physics Exchange Program for Studies of Hadron Structure and QCD) by the US Department of Energy under grant DE-SC0006758.

Appendix A BARYONS AS QUARK - DIQUARK BOUND STATES

The quark - diquark model, based on the static approximation to the Faddeev equation, for octet baryons in the limit of isospin symmetry ( $M_{u}=M_{d}$ ) has been described in Refs. Carrillo-Serrano et al. (2014, 2016). As explained in the main text, in our present work we still assume isospin symmetry in the vacuum, and therefore equal current quark masses ( $m_{u}=m_{d}\equiv m$ ) and equal constituent quark masses in the vacuum ( $M_{u0}=M_{d0}\equiv M_{0}$ ). However, a consistent description of isospin asymmetric systems, like neutron star matter, in the framework of an effective quark theory requires to consider the spontaneous breaking of isospin symmetry due to the presence of the medium, i.e., $M_{u}\neq M_{d}$ . In this appendix, we therefore briefly explain the main points of our model for the octet baryons, treating the masses $M_{u},M_{d},M_{s}$ as independent quantities.

The chiral invariant interaction Lagrangian in the $qq$ channel is given by Vogl and Weise (1991)

$\displaystyle{\cal L}_{I}^{(qq)}$	$\displaystyle=G_{S}\left[\left(\bar{q}\gamma_{5}\,C\,\lambda_{a}\,\lambda^{(C)}_{A}\,\bar{q}^{T}\right)\,\left({q}^{T}C^{-1}\gamma_{5}\,\lambda_{a}\,\lambda^{(C)}_{A}\,q\right)\right.$
	$\displaystyle\left.-\left(\bar{q}\,C\,\lambda_{a}\,\lambda^{(C)}_{A}\,\bar{q}^{T}\right)\,\left({q}^{T}C^{-1}\,\lambda_{a}\,\lambda^{(C)}_{A}\,q\right)\right]$
	$\displaystyle+G_{A}\left[\left(\bar{q}\gamma_{\mu}\,C\,\lambda_{s}\,\lambda^{(C)}_{A}\,\bar{q}^{T}\right)\,\left({q}^{T}C^{-1}\gamma_{\mu}\,\lambda_{s}\,\lambda^{(C)}_{A}\,q\right)\right.$
	$\displaystyle\left.+\left(\bar{q}\gamma_{\mu}\gamma_{5}\,C\,\lambda_{a}\,\lambda^{(C)}_{A}\,\bar{q}^{T}\right)\,\left({q}^{T}C^{-1}\,\gamma^{\mu}\gamma_{5}\,\lambda_{a}\,\lambda^{(C)}_{A}\,q\right)\right]\,.$	(49)

Here $\lambda_{a}$ ( $a=2,5,7)$ are the antisymmetric Gell-Mann flavor matrices, $\lambda_{s}$ ( $s=0,1,3,4,6,8$ ) are the symmetric ones, and the antisymmetric Gell-Mann color matrices $\lambda^{(C)}_{A}$ ( $A=2,5,7$ ) project to color $\bar{3}$ diquark states. (There are also interaction terms in the color $6$ diquark channels, which are not shown here because they do not contribute to colorless baryon states.) The charge conjugation Dirac matrix is $C=i\gamma_{2}\gamma_{0}$ . The first line in (49) is the interaction in the scalar diquark ( $0^{+}$ ) channel, the second line shows the pseudoscalar diquark ( $0^{-}$ ) channel, the third line the axial vector diquark ( $1^{+}$ ) channel, and the fourth line the vector diquark ( $1^{-}$ ) channel. Following previous works Carrillo-Serrano et al. (2014, 2016), we will include only the scalar and the axial vector diquark channels, which are expected to be dominant from the nonrelativistic analogy.

By simple manipulations in flavor space, we can identically rewrite the 2 terms relevant for our calculation as follows:

	$\displaystyle{\cal L}^{(qq)}$	$\displaystyle=G_{S}\left(\bar{q}\gamma_{5}\,C\,t_{a}\,\lambda^{(C)}_{A}\,\bar{q}^{T}\right)\,\left({q}^{T}C^{-1}\gamma_{5}\,t_{a}^{\dagger}\,\lambda^{(C)}_{A}\,q\right)$
		$\displaystyle+G_{A}\left(\bar{q}\gamma_{\mu}\,C\,t_{s}\,\lambda^{(C)}_{A}\,\bar{q}^{T}\right)\,\left({q}^{T}C^{-1}\gamma_{\mu}\,t_{s}^{\dagger}\,\lambda^{(C)}_{A}\,q\right)\,,$		(50)

where we introduced the three antisymmetric and anti-Hermitian $3\times 3$ flavor matrices $t_{a}\equiv\left(t_{[ud]},t_{[us]},t_{[ds]}\right)$ , and the six symmetric and Hermitian $3\times 3$ flavor matrices $t_{s}\equiv\left(t_{\{ud\}},t_{\{us\}},t_{\{ds\}},t_{\{uu\}},t_{\{dd\}},t_{\{ss\}}\right)$ . For example, $t_{[ud]}$ is given by $\displaystyle{t_{[ud]}=\left(\begin{array}[]{ccc}0&1&0\\ -1&0&0\\ 0&0&0\end{array}\right)}$ , and corresponds to the antisymmetric flavor combination expressed by $[ud]$ . The matrices $t_{[us]}$ , $t_{[ds]}$ are defined similarly, and correspond to the antisymmetric flavor combinations $[us]$ and $[ds]$ . The symmetric matrix $t_{\{ud\}}$ has the same structure as $t_{[ud]}$ , but with the $-1$ replaced by $+1$ in the (2,1) component, and corresponds to the symmetric flavor combination expressed by $\{ud\}$ . The matrices $t_{\{us\}}$ and $t_{\{ds\}}$ are defined in a similar way for the symmetric flavor combinations $\{us\}$ and $\{ds\}$ . Finally, the matrices $t_{\{uu\}}$ , $t_{\{dd\}}$ and $t_{\{ss\}}$ have a $\sqrt{2}$ as the (1,1) component, the (2,2) component, and the (3,3) component, respectively, with all other components equal to zero.

We express the Faddeev vertex functions for a given baryon by $X^{a}_{i}$ , where $a$ denotes the diquark channels explained above ( $a=[ud],\dots,\{ss\}$ ), and $i$ is the flavor of the third quark. For example, for the proton the diquark-quark channels are labeled by $[ud]u$ , $\{ud\}u$ , and $\{uu\}d$ . The Faddeev equations for the vertex functions $X^{a}_{i}(p,q)$ , describing a baryon of momentum $p$ as a bound state of a quark (momentum $q$ ) and a diquark (momentum $p-q$ ), are

\displaystyle X^{a}_{i}(p,q)=\int\frac{{\rm d}^{4}k}{(2\pi)^{4}}\,Z^{ab}_{ij}\,S_{j}(k)\,\tau^{bc}_{(ki)}(p-k)\,X^{c}_{j}(p,k)\,,

(51)

which is shown graphically in Fig. 10. Here the quark exchange kernel is given by

\displaystyle Z^{ab}_{ij}=-3\Lambda^{b}\,\left(t^{b}\,S(k+q-p)t^{a\dagger}\right)_{ij}\,\Lambda^{a}\,,

(52)

where the factor $-3$ comes from projection to color singlet states, and we used the identity $CS^{T}(k)C^{-1}=S(-k)$ to process the charge conjugation matrices. The Dirac matrices $\Lambda$ are given by $\Lambda^{a}=\gamma_{5}$ for the scalar diquark channels (flavor index $a=[ud],[us],[ds]$ ) and $\Lambda^{a}=\gamma^{\mu}\gamma_{5}$ for the axial vector diquark channels (flavor index $a=\{ud\},\{us\},\{ds\},\{uu\},\{dd\},\{ss\}$ ). The quantities $\tau^{bc}_{(ki)}$ in Eq. (51) are diagonal in the diquark flavor indices, $\tau^{bc}_{(ki)}=\delta_{bc}\tau^{b}_{(ki)}$ , where $\tau^{b}_{(ki)}$ is the reduced $t$ -matrix in the diquark channel $b$ with interacting quark flavors $k,i$ . Therefore, $\tau^{b}_{(ki)}\equiv\tau_{[ki]}$ in the scalar diquark channels ( $[ki]=[ud],[us],[ds]$ ), and $\tau^{b}_{(ki)}\equiv\tau^{\mu\nu}_{\{ki\}}$ in the axial vector diquark channels ( $\{ki\}=\{ud\},\{us\},\{ds\},\{uu\},\{dd\},\{ss\}$ ). The explicit forms of the reduced diquark $t$ -matrices $\tau_{[ki]}$ and $\tau^{\mu\nu}_{\{ki\}}$ are given in Ref. Carrillo-Serrano et al. (2014, 2016). In the last factor $\left(\dots\right)$ of (52), the quark propagator is considered as a $3\times 3$ diagonal matrix with diagonal elements $S_{k}=(S_{u},S_{d},S_{s})$ .

The formalism described so far is the Faddeev framework in the NJL model, where the only assumptions are the ladder approximation for the 2-body $t$ -matrices and the restriction to the scalar and axial vector diquark channels. The quark-diquark model used in the calculations of the main text is obtained by the replacement $S\rightarrow-\frac{1}{M}$ in the quark exchange kernel (52) for each quark flavor $k$ , i.e., this approximation neglects the momentum dependence of the quark exchange kernel and is called the “static approximation” of the Faddeev kernel Hellstern et al. (1997). In this approximation, the vertex functions $X^{a}_{i}$ of (51) depend only on the total momentum $p$ , and the integral in Eq. (51) extends only over the product $S_{j}(k)\,\tau^{bc}_{(ki)}(p-k)$ , which we regularize according to the proper time scheme (see App. C), which avoids unphysical thresholds for the decay into quarks. The Dirac structure of the vertex function can also be determined analytically. Let us again take the proton as an example. Arranging the three interacting channels mentioned above into a vector, the vertex function can be expressed in the form ¹³¹³13The coupling of the time component of the axial vector diquark and the quark to the total spin $(\frac{1}{2},S_{N})$ gives rise to the structure $\frac{p^{\mu}}{M_{N}}\gamma_{5}u(p,S_{N})$ , and the coupling of the 3-vector components of the axial vector diquark to the quark gives the structure $\displaystyle\sum_{\lambda,s}\left(1\frac{1}{2},\lambda s|\frac{1}{2}S_{N}\right)\varepsilon^{\mu}(p,\lambda)\,u(p,s)=\frac{-1}{\sqrt{3}}\left(\frac{p^{\mu}}{M_{N}}+\gamma^{\mu}\right)\gamma_{5}\,u(p,S_{N})\,,$ where $\varepsilon^{\mu}(p,\lambda)$ is the Lorentz 4-vector for spin 1 with mass $M_{N}$ , and $u(p,S_{N})$ is the Dirac spinor with mass $M_{N}$ .

\displaystyle|p\rangle=\left(\begin{array}[]{c}X^{[ud]}_{u}(p)\\ X^{\{ud\}}_{u}(p)\\ X^{\{uu\}}_{d}(p)\end{array}\right)\equiv\left(\begin{array}[]{c}\alpha_{1}\,[ud]u\\ \left(\alpha_{2}\frac{p^{\mu}}{M_{p}}+\alpha_{3}\gamma^{\mu}\right)\gamma_{5}\,\{ud\}u\\ \left(\alpha_{4}\frac{p^{\mu}}{M_{p}}+\alpha_{5}\gamma^{\mu}\right)\gamma_{5}\,\{uu\}d\end{array}\right)\,u_{p}(p)\,,

(59)

where $u_{p}(p)$ is the Dirac spinor with the mass of the proton ( $M_{p}$ ). Inserting (59) into Eq. (51) then gives homogeneous equations for the coefficients $\alpha_{i}$ , and the characteristic equation gives the proton mass $M_{p}=M_{p}(M_{u},M_{d})$ .

The vertex functions of the other members of the octet with two identical quark flavors are similar to (59), with obvious replacements of quark flavors. For the $\Sigma^{0}$ , the flavor structure can be obtained by acting with the isospin lowering operator ( $T_{-}$ ) on $|\Sigma^{+}\rangle$ , which generates

\displaystyle|\Sigma^{0}\rangle=\left(\begin{array}[]{c}X^{[us]}_{d}(p)\\ X^{[ds]}_{u}(p)\\ X^{\{us\}}_{d}(p)\\ X^{\{ds\}}_{u}(p)\\ X^{\{ud\}}_{s}(p)\end{array}\right)\equiv\left(\begin{array}[]{c}\alpha_{1}\,[us]d\\ \alpha_{2}\,[ds]u\\ \left(\alpha_{3}\frac{p^{\mu}}{M_{\Sigma^{0}}}+\alpha_{4}\gamma^{\mu}\right)\gamma_{5}\,\{us\}d\\ \left(\alpha_{5}\frac{p^{\mu}}{M_{\Sigma^{0}}}+\alpha_{6}\gamma^{\mu}\right)\gamma_{5}\,\{ds\}u\\ \left(\alpha_{7}\frac{p^{\mu}}{M_{\Sigma^{0}}}+\alpha_{8}\gamma^{\mu}\right)\gamma_{5}\,\{ud\}s\end{array}\right)\,u_{\Sigma^{0}}(p)\,,

(70)

where $u_{\Sigma^{0}}(p)$ is the Dirac spinor with the mass $M_{\Sigma^{0}}$ . Note that there is no component with the flavor structure $[ud]s$ in the $\Sigma^{0}$ , and, of course, also no components where the 2 light quarks form a scalar diquark in $\Sigma^{\pm}$ because those vanish identically ( $[uu]=[dd]=0$ ).

For the $\Lambda$ , we first construct a state $U_{+}|\Xi^{0}\rangle$ , where the raising $U$ -spin operator converts an $s$ -quark into a $d$ -quark, and orthogonalize this state to $|\Sigma^{0}\rangle$ . This gives

\displaystyle|\Lambda\rangle=\left(\begin{array}[]{c}X^{[ud]}_{s}(p)\\ X^{[us]}_{d}(p)\\ X^{[ds]}_{u}(p)\\ X^{\{us\}}_{d}(p)\\ X^{\{ds\}}_{u}(p)\end{array}\right)\equiv\left(\begin{array}[]{c}\alpha_{1}\,[ud]s\\ \alpha_{2}\,[us]d\\ \alpha_{3}\,[ds]u\\ \left(\alpha_{4}\frac{p^{\mu}}{M_{\Lambda}}+\alpha_{5}\gamma^{\mu}\right)\gamma_{5}\,\{us\}d\\ \left(\alpha_{6}\frac{p^{\mu}}{M_{\Lambda}}+\alpha_{7}\gamma^{\mu}\right)\gamma_{5}\,\{ds\}u\end{array}\right)\,u_{\Lambda}(p)\,,

(81)

where $u_{\Lambda}(p)$ is the Dirac spinor with the mass $M_{\Lambda}$ . Note that there is no component with the flavor structure $\{ud\}s$ in the $\Lambda$ .

In the calculations of the main text we only need the masses of the octet baryons as functions of the constituent quark masses. In isospin asymmetric systems like neutron star matter, the isospin symmetry is obviously broken, but the charge symmetry is intact if we simultaneously reverse the signs of the isospin $z$ -components of the baryons and the constituent quarks in the baryon. We therefore have the following five independent functions of $M_{u},M_{d}$ (omitting the obvious dependence on $M_{s}$ for the hyperons to simplify the notations):

	$\displaystyle M_{p}$	$\displaystyle=M_{p}(M_{u},M_{d})\,,\,\,\,\,\,\,\,M_{\Sigma^{+}}=M_{\Sigma^{+}}(M_{u})\,,\,\,\,\,\,\,M_{\Xi^{+}}=M_{\Xi^{+}}(M_{u})\,,$
	$\displaystyle M_{\Sigma^{0}}$	$\displaystyle=M_{\Sigma^{0}}(M_{u},M_{d})=M_{\Sigma^{0}}(M_{d},M_{u})\,,$
	$\displaystyle M_{\Lambda}$	$\displaystyle=M_{\Lambda}(M_{u},M_{d})=M_{\Lambda}(M_{d},M_{u})\,.$

The masses of the remaining baryons can then be expressed by

	$\displaystyle M_{n}$	$\displaystyle=M_{p}(M_{d},M_{u})\,,$	$\displaystyle M_{\Sigma^{-}}$	$\displaystyle=M_{\Sigma^{+}}(M_{d})\,,$
	$\displaystyle M_{\Xi^{0}}$	$\displaystyle=M_{\Xi^{+}}(M_{d})\,.$				(82)

The vertex functions and masses of the decuplet baryons are calculated in a similar way. The calculation is simplified by the fact that here only the axial vector diquark channels (symmetric combinations of quark flavors $\{q_{1}q_{2}\}$ ) contribute, which leaves only one possible Dirac-Lorentz structure for all components of the baryon vertex, namely the Rarita-Schwinger spinor $u^{\mu}(p,S_{b})$ .

In the present work, we determined the coupling constants $G_{S}$ , $G_{A}$ of the Lagrangian (49) so as to reproduce the observed masses of the nucleon ( $M_{N}=0.94$ GeV) and the Delta baryon ( $M_{\Delta}=1.232$ GeV) in the vacuum ( $M_{u}=M_{d}=0.4$ GeV). We also determined our vacuum value of the strange quark mass so as to reproduce the observed mass of the $\Omega$ baryon ( $M_{\Omega}=1.67$ GeV). In this way we obtain

\displaystyle G_{S}

\displaystyle=8.76\,{\rm GeV}^{-2}\,,

\displaystyle G_{A}

\displaystyle=7.36\,{\rm GeV}^{-2}\,,

\displaystyle M_{s}

\displaystyle=0.562\,{\rm GeV}\,.

(83)

The resulting masses of octet baryons are given in Tab. 2 of the main text. The masses of the decuplet baryons - except for the $\Delta$ and the $\Omega$ which were fitted - are also well reproduced in this calculation; we obtain $M_{\Sigma^{*}}=1.38$ GeV, and $M_{\Xi^{*}}=1.53$ GeV. The mass of the kaon, however, is underestimated ( $0.43$ GeV for the case of 4-fermi couplings); in order to reproduce its observed mass we would need a larger value of $M_{s}$ . Because the focus of our present work is on the baryons, we made no attempt to reproduce the meson masses well.

The masses of the diquarks (poles of the quantities $\tau$ in Eq. (51)) are $M^{[\ell\ell^{\prime}]}=0.768$ GeV and $M^{[\ell s]}=0.902$ GeV for the scalar diquarks with $\ell,\ell^{\prime}=u,d$ , and $M^{\{\ell\ell^{\prime}\}}=0.929$ GeV, $M^{\{\ell s\}}=1.04$ GeV, $M^{\{ss\}}=1.15$ GeV for the axial vector diquarks.

We finally mention that the values of $G_{S}$ and $G_{A}$ given in (83) are different from those used in a previous work on the flavor $SU(2)$ case Tanimoto et al. (2020). There $G_{S}$ and $G_{A}$ were fitted to the mass of the free nucleon and its axial vector coupling constant ( $g_{A}=1.26$ ). The values of $G_{S}$ ( $G_{A}$ ) obtained in that way were larger (smaller) than the values given in (83), which indicates that the dominance of the scalar diquark channel, which increases in-medium because of the decreasing scalar diquark mass, was more pronounced in Ref. Tanimoto et al. (2020) than in the present work. This stronger attraction in-medium, however, was eventually canceled by a stronger repulsion in the vector-isovector $q\bar{q}$ channel, because in the flavor $SU(2)$ case it was possible to reproduce the symmetry energy in the mean field approximation without violating the chiral symmetry of the interaction Lagrangian. As a result, the pressure in neutron star matter and the star masses calculated in Ref. Tanimoto et al. (2020) were almost identical to the results shown by the dashed lines in Fig. 8 of the main text.

In our present work, we found it more essential to reproduce the $N-\Delta$ mass difference, because one of our motivations was to see how the $\Sigma-\Lambda$ mass difference evolves with density if the spin dependent diquark correlations are constrained to the vacuum value of the $N-\Delta$ mass difference from the outset. The value of $g_{A}$ , obtained with our present coupling constants (83), is larger than the observed value by about 10 $\%$ . In future investigations, we wish to see whether the inclusion of additional diquark channels allows one to find a set of coupling constants which reproduces the nucleon mass, the Delta baryon mass, and $g_{A}$ simultaneously.

Appendix B Meson exchange in symmetric nuclear matter

As mentioned in Sec. III.1, in order to express the effective baryon-nucleon interactions (17) and (18) in terms of meson exchange processes of the type shown in Fig. 11, one should multiply the numerator and denominator functions in the first two lines of those expressions by the quark-meson couplings

\displaystyle g_{\sigma}^{(q)2}=g_{\delta}^{(q)2}=\frac{-1}{\Pi^{\prime}_{s}(q^{2}=0)}\,,

(84)

and similarly in the third lines by

\displaystyle g_{\omega}^{(q)2}=g_{\rho}^{(q)2}=\frac{-1}{\Pi^{\prime}_{v}(q^{2}=0)}\,.

(85)

Here the $\overline{q}q$ bubble graphs in the scalar and vector channels are given by

$\displaystyle\Pi_{s}(q^{2})$	$\displaystyle=12i\int\frac{{\rm d}^{4}k}{(2\pi)^{4}}$
	$\displaystyle\times\left[\frac{-2}{k^{2}-M^{2}}+\left(q^{2}-4M^{2}\right)\int_{0}^{1}\,{\rm d}x\frac{1}{\left(k^{2}-M^{2}+q^{2}x(1-x)\right)^{2}}\right]\,,$	(86)
$\displaystyle\Pi_{v}(q^{2})$	$\displaystyle=48iq^{2}\,\int\frac{{\rm d}^{4}k}{(2\pi)^{4}}\,\int_{0}^{1}\,{\rm d}x\frac{x(1-x)}{\left(k^{2}-M^{2}+q^{2}x(1-x)\right)^{2}}\,,$	(87)

and the primes in (84) and (85) mean differentiation w.r.t. $q^{2}$ .

Table 6: Effective coupling constants and masses of

\sigma

\omega

\delta

and

\rho

mesons for four values of the baryon density in symmetric nuclear matter. Coupling constants are dimensionless, and masses are given in units of GeV. For definitions, see Eqs. (88), (89), and text.

density	$g_{\sigma}^{(q)}$	$g_{\sigma}^{(N)}$	$g_{\sigma}^{(\Lambda)}$	$g_{\sigma}^{(\Sigma)}$	$g_{\sigma}^{(\Xi)}$	$M_{\sigma}$
0	6.63	18.21	12.13	10.28	5.70	1.25
0.15	4.85	10.06	7.05	5.59	3.15	0.96
0.3	4.20	6.25	4.80	3.59	2.06	0.97
0.5	3.85	4.01	3.45	2.44	1.43	1.05
density	$g_{\omega}^{(q)}$	$g_{\omega}^{(N)}$	$g_{\omega}^{(\Lambda)}$	$g_{\omega}^{(\Sigma)}$	$g_{\omega}^{(\Xi)}$	$M_{\omega}$
0	5.27	15.80	10.54	10.54	5.27	1.52
0.15	4.51	13.53	9.02	9.02	4.51	1.30
0.3	4.18	12.53	8.35	8.35	4.18	1.20
0.5	3.99	11.96	7.97	7.97	3.99	1.15
density	$g_{\delta}^{(q)}$	$g_{\delta}^{(p)}$	$g_{\delta}^{(\Lambda)}$	$g_{\delta}^{(\Sigma^{+})}$	$g_{\delta}^{(\Xi^{0})}$	$M_{\delta}$
0	6.63	4.64	0	10.28	5.70	1.25
0.15	4.85	2.38	0	5.59	3.15	0.99
0.3	4.20	1.38	0	3.59	2.06	1.00
0.5	3.85	0.81	0	2.44	1.43	1.10
density	$g_{\rho}^{(q)}$	$g_{\rho}^{(p)}$	$g_{\rho}^{(\Lambda)}$	$g_{\rho}^{(\Sigma^{+})}$	$g_{\rho}^{(\Xi^{0})}$	$M_{\rho}$
0	5.27	5.27	0	10.54	5.27	1.52
0.15	4.51	4.51	0	9.02	4.51	1.30
0.3	4.18	4.18	0	8.35	4.18	1.20
0.5	3.99	3.99	0	7.97	3.99	1.15

For simplicity we consider only the $\ell=0$ terms in (17) and (18). They can be expressed in the following form

	$\displaystyle f_{0,bN}$	$\displaystyle=-\frac{M_{b}}{E_{b}}\frac{M_{N}}{E_{N}}\,\frac{g_{\sigma}^{(b)}g_{\sigma}^{(N)}}{M_{\sigma}^{2}}+\frac{g_{\omega}^{(b)}g_{\omega}^{(N)}}{M_{\omega}^{2}}\,,$		(88)
	$\displaystyle f^{\prime}_{0,bN}$	$\displaystyle=-\frac{M_{b}}{E_{b}}\frac{M_{N}}{E_{N}}\,\frac{g_{\delta}^{(b)}g_{\delta}^{(p)}}{M_{\delta}^{2}}+\frac{g_{\rho}^{(b)}g_{\rho}^{(p)}}{M_{\rho}^{2}}\,.$		(89)

Here all meson-baryon coupling constants and meson masses are defined at zero momentum, and are different from the values at the meson poles. The resulting values for the effective coupling constants and masses are summarized for three values of the baryon density in Tab. 6. In relation to our discussions in Sec. III.3, we note that $g_{\sigma}^{(\Lambda)}>g_{\sigma}^{(\Sigma)}$ , which reflects the different internal quark-diquark structure of the $\Lambda$ and the $\Sigma$ baryons.

We finally add a few comments on the definition of effective coupling constants and meson masses used here: First, the multiplication of the density dependent scaling factors (84) and (85) to the numerators and denominators of (17) and (18) obscures the simplicity of those basic expressions, and for better orientation the values listed in Tab. 4 of the main text is more useful. Nevertheless, it is necessary for a proper definition of coupling constants and meson masses at zero momentum of the mesons. For the coupling constants, this is immediately clear from Fig. 11. For the meson masses, consider for example the case of the $\sigma$ meson. The reduced $t$ -matrix in the $0^{+}$ $\overline{q}q$ channel is given by

\displaystyle\tau_{\sigma}(q^{2})=\frac{-2G_{\pi}}{1+2G_{\pi}\Pi_{s}(q^{2})+2G_{\pi}\,\delta M_{\sigma}^{2}}\,,

(90)

where the nucleon loop contributions, approximated by their forms at $q=0$ , are denoted by $\delta M_{\sigma}^{2}$ . Expanding (90) around $q^{2}=0$ gives the approximate Yukawa-like form

\displaystyle\tau_{\sigma}(q^{2})=\frac{g_{\sigma}^{(q)2}}{q^{2}-M_{\sigma}^{2}}\,,

where $g_{\sigma}^{(q)}$ is defined by (84), and

\displaystyle M_{\sigma}^{2}=g_{\sigma}^{(q)2}\left(\frac{1}{2G_{\pi}}+\Pi_{s}(q^{2}=0)+\delta M_{\sigma}^{2}\right)\,.

(91)

The terms $\left(\dots\right)$ in (91) agree with the denominator in the second line of Eq. (17) because of the relation $\Pi_{s}(q^{2}=0)=2g(M)$ , where $g(M)$ is given by Eq. (19).

Appendix C Regularization method

To evaluate 4-dimensional integrals, we introduce Feynman parameters and perform shifts of the loop momentum so that the integrand depends only on $k^{2}$ , where $k$ is the loop momentum, besides other fixed variables. We then perform a Wick rotation and use 4-dimensional spherical polar coordinates to obtain

\displaystyle\int\mathop{}\!\mathrm{d}^{4}k\,f(k^{2})=2\pi^{2}i\int^{\infty}_{0}\mathop{}\!\mathrm{d}k_{E}\,k_{E}^{3}\,f(-k_{E}^{2})\,,

where $k_{E}=\sqrt{k_{0}^{2}+\bm{k}^{2}}$ is the Euclidean length. Next, we consider the following identities:

	$\displaystyle\ln\frac{D}{D_{0}}$	$\displaystyle=-\int_{0}^{\infty}\frac{{\rm d}\tau}{\tau}\left(e^{-\tau D}-e^{-\tau D_{0}}\right)\,,$		(92)
	$\displaystyle\frac{1}{D^{n}}$	$\displaystyle=\frac{1}{(n-1)!}\int^{\infty}_{0}\,{\rm d}\tau\,\tau^{n-1}e^{-\tau D}\,\,\,\,\,\,\,\,\,\,\,\,(n\geq 0)\,,$		(93)

where $D$ is a function of $k_{E}^{2}$ and other fixed variables. In the proper time regularization scheme, the infrared cutoff ( $\Lambda_{\rm IR}$ ) is introduced by replacing the upper integration limits in (92), (93) by $1/\Lambda_{\rm IR}^{2}$ , and the ultraviolet cutoff ( $\Lambda_{\rm UV}$ ) by replacing the lower integration limits by $1/\Lambda_{\rm UV}^{2}$ . After these replacements, one performs the integration over $k_{E}$ . The ultraviolet cutoff makes the integrals finite, while the infrared cutoff eliminates unphysical thresholds (imaginary parts) for the decay of hadrons into quarks, thus simulating the role of confinement.

Appendix D Sizes of quark cores in the nuclear medium

The rms radius of the baryon density distribution of the quark core of the nucleons in the medium is related to the isoscalar combination of the corresponding electric charge radii of protons and neutrons by

\displaystyle r_{N}(\rho_{B})=\sqrt{\langle r_{Ep}^{2}\rangle(\rho_{B})+\langle r_{En}^{2}\rangle(\rho_{B})}\,.

(94)

In the language of Feynman diagrams used in Ref.Cloët et al. (2014), the corresponding isoscalar baryon form factor is obtained by the operator insertion $\frac{1}{3}\gamma^{0}$ on each quark line. For free nucleons (zero density) the result of the NJL model calculations of Ref. Cloët et al. (2014), using the same parameters as in the present paper, is $r_{N}(0)=0.475$ fm. Note that this is the value for the quark core without meson cloud corrections, obtained by replacing the dressed quark form factors in Sect. VI of Ref. Cloët et al. (2014) by their bare values ( $F_{1U}=\frac{2}{3},F_{1D}=-\frac{1}{3},F_{2U}=F_{2D}=0$ ). The pion cloud contributions to the isoscalar quantity $r_{N}$ are very small. A simple estimate of $\omega$ meson cloud effects, using our present value of $G_{v}$ , gives only a small correction, but a more realistic treatment, following the lines of the vector meson dominance model with the observed $\omega$ meson pole, increases the isoscalar baryon radius to $0.78$ fm (see Table VI of Ref. Cloët et al. (2014)), which is close to the experimental value. As mentioned in the main text, however, the quantity which seems more relevant for the role of the Pauli principle is the baryon radius of the quark core without meson cloud effects. This is simply because the mesons are bosons, and the overlap of the meson clouds just corresponds to the meson exchange interactions. Therefore the results shown in Sec. III.3.5 of the main text refer to this quantity.

Based on naive geometric intuition, the volume fraction occupied by the quark cores can be defined as

\displaystyle f(\rho_{B})=\rho_{B}\,v_{N}(\rho_{B})\,,

(95)

where $v_{N}(\rho_{B})=\frac{4\pi}{3}r_{N}^{3}(\rho_{B})$ . Estimates based on this expression are also given in the main text.

References

Haidenbauer and Meissner (2005) J. Haidenbauer and U.-G. Meissner, Phys. Rev. C 72, 044005 (2005), arXiv:nucl-th/0506019 .
Rijken and Yamamoto (2006) T. A. Rijken and Y. Yamamoto, Phys. Rev. C 73, 044008 (2006), arXiv:nucl-th/0603042 .
Bogner et al. (2010) S. K. Bogner, R. J. Furnstahl, and A. Schwenk, Prog. Part. Nucl. Phys. 65, 94 (2010), arXiv:0912.3688 [nucl-th] .
Petschauer et al. (2020) S. Petschauer, J. Haidenbauer, N. Kaiser, U.-G. Meißner, and W. Weise, Front. in Phys. 8, 12 (2020), arXiv:2002.00424 [nucl-th] .
Lonardoni et al. (2014) D. Lonardoni, F. Pederiva, and S. Gandolfi, Phys. Rev. C 89, 014314 (2014), arXiv:1312.3844 [nucl-th] .
Faessler et al. (1982) A. Faessler, F. Fernandez, G. Lubeck, and K. Shimizu, Phys. Lett. B 112, 201 (1982).
Oka et al. (1987) M. Oka, K. Shimizu, and K. Yazaki, Nucl. Phys. A 464, 700 (1987).
Fujiwara et al. (2007) Y. Fujiwara, Y. Suzuki, and C. Nakamoto, Prog. Part. Nucl. Phys. 58, 439 (2007), arXiv:nucl-th/0607013 .
Oka (2023) M. Oka, (2023), arXiv:2301.06026 [nucl-th] .
Inoue et al. (2010) T. Inoue, N. Ishii, S. Aoki, T. Doi, T. Hatsuda, Y. Ikeda, K. Murano, H. Nemura, and K. Sasaki (HAL QCD), Prog. Theor. Phys. 124, 591 (2010), arXiv:1007.3559 [hep-lat] .
Hyodo and Niiyama (2021) T. Hyodo and M. Niiyama, Prog. Part. Nucl. Phys. 120, 103868 (2021), arXiv:2010.07592 [hep-ph] .
Gibson et al. (2010) B. F. Gibson, K. Imai, T. Motoba, T. Nagae, and A. Ohnishi, eds., Proceedings, 10th International Conference on Hypernuclear and strange particle physics (HYP 2009): Tokai, Japan, September 14-18, 2009, Vol. 835 (2010).
Gal et al. (2016) A. Gal, E. V. Hungerford, and D. J. Millener, Rev. Mod. Phys. 88, 035004 (2016), arXiv:1605.00557 [nucl-th] .
Hiyama and Nakazawa (2018) E. Hiyama and K. Nakazawa, Ann. Rev. Nucl. Part. Sci. 68, 131 (2018).
Tamura (2022) H. Tamura, EPJ Web Conf. 271, 12001 (2022).
Miwa et al. (2022) K. Miwa et al., EPJ Web Conf. 271, 04001 (2022).
Nanamura et al. (2022) T. Nanamura et al. (J-PARC E40), PTEP 2022, 093D01 (2022), arXiv:2203.08393 [nucl-ex] .
Demorest et al. (2010) P. Demorest, T. Pennucci, S. Ransom, M. Roberts, and J. Hessels, Nature 467, 1081 (2010), arXiv:1010.5788 [astro-ph.HE] .
Antoniadis et al. (2013) J. Antoniadis et al., Science 340, 6131 (2013), arXiv:1304.6875 [astro-ph.HE] .
Riley et al. (2021) T. E. Riley et al., Astrophys. J. Lett. 918, L27 (2021), arXiv:2105.06980 [astro-ph.HE] .
Fonseca et al. (2021) E. Fonseca et al., Astrophys. J. Lett. 915, L12 (2021), arXiv:2104.00880 [astro-ph.HE] .
Glendenning (1997) N. K. Glendenning, Compact stars: Nuclear physics, particle physics, and general relativity (1997).
Bombaci (2017) I. Bombaci, JPS Conf. Proc. 17, 101002 (2017), arXiv:1601.05339 [nucl-th] .
Chatterjee and Vidaña (2016) D. Chatterjee and I. Vidaña, Eur. Phys. J. A 52, 29 (2016), arXiv:1510.06306 [nucl-th] .
Burgio et al. (2021) G. F. Burgio, H. J. Schulze, I. Vidana, and J. B. Wei, Prog. Part. Nucl. Phys. 120, 103879 (2021), arXiv:2105.03747 [nucl-th] .
Weissenborn et al. (2012) S. Weissenborn, D. Chatterjee, and J. Schaffner-Bielich, Phys. Rev. C 85, 065802 (2012), [Erratum: Phys.Rev.C 90, 019904 (2014)], arXiv:1112.0234 [astro-ph.HE] .
Spinella and Weber (2019) W. M. Spinella and F. Weber, Astron. Nachr. 340, 145 (2019), arXiv:1812.03600 [nucl-th] .
Yamamoto et al. (2013) Y. Yamamoto, T. Furumoto, N. Yasutake, and T. A. Rijken, Phys. Rev. C 88, 022801 (2013), arXiv:1308.2130 [nucl-th] .
Haidenbauer et al. (2017) J. Haidenbauer, U. G. Meißner, N. Kaiser, and W. Weise, Eur. Phys. J. A 53, 121 (2017), arXiv:1612.03758 [nucl-th] .
Kohno (2018) M. Kohno, Phys. Rev. C 97, 035206 (2018), arXiv:1802.05388 [nucl-th] .
Logoteta et al. (2019) D. Logoteta, I. Vidana, and I. Bombaci, Eur. Phys. J. A 55, 207 (2019), arXiv:1906.11722 [nucl-th] .
Gerstung et al. (2020) D. Gerstung, N. Kaiser, and W. Weise, Eur. Phys. J. A 56, 175 (2020), arXiv:2001.10563 [nucl-th] .
Contrera et al. (2022) G. A. Contrera, D. Blaschke, J. P. Carlomagno, A. G. Grunfeld, and S. Liebing, Phys. Rev. C 105, 045808 (2022), arXiv:2201.00477 [nucl-th] .
Buballa (2005) M. Buballa, Phys. Rept. 407, 205 (2005), arXiv:hep-ph/0402234 [hep-ph] .
Alford et al. (2008) M. G. Alford, A. Schmitt, K. Rajagopal, and T. Schäfer, Rev. Mod. Phys. 80, 1455 (2008), arXiv:0709.4635 [hep-ph] .
Guichon (1988) P. A. M. Guichon, Phys. Lett. B 200, 235 (1988).
Chodos et al. (1974) A. Chodos, R. L. Jaffe, K. Johnson, C. B. Thorn, and V. F. Weisskopf, Phys. Rev. D 9, 3471 (1974).
Nambu and Jona-Lasinio (1961a) Y. Nambu and G. Jona-Lasinio, Phys. Rev. 122, 345 (1961a), [,127(1961)].
Nambu and Jona-Lasinio (1961b) Y. Nambu and G. Jona-Lasinio, Phys. Rev. 124, 246 (1961b), [,141(1961)].
Vogl and Weise (1991) U. Vogl and W. Weise, Prog. Part. Nucl. Phys. 27, 195 (1991).
Hatsuda and Kunihiro (1994) T. Hatsuda and T. Kunihiro, Phys. Rept. 247, 221 (1994), arXiv:hep-ph/9401310 [hep-ph] .
Ishii et al. (1995) N. Ishii, W. Bentz, and K. Yazaki, Nucl. Phys. A 587, 617 (1995).
Bentz and Thomas (2001) W. Bentz and A. W. Thomas, Nucl. Phys. A696, 138 (2001), arXiv:nucl-th/0105022 [nucl-th] .
Geesaman et al. (1995) D. F. Geesaman, K. Saito, and A. W. Thomas, Ann. Rev. Nucl. Part. Sci. 45, 337 (1995).
Lu et al. (1998) D.-H. Lu, A. W. Thomas, K. Tsushima, A. G. Williams, and K. Saito, Phys. Lett. B 417, 217 (1998), arXiv:nucl-th/9706043 .
Stone et al. (2016) J. R. Stone, P. A. M. Guichon, P. G. Reinhard, and A. W. Thomas, Phys. Rev. Lett. 116, 092501 (2016), arXiv:1601.08131 [nucl-th] .
Cloët et al. (2005) I. C. Cloët, W. Bentz, and A. W. Thomas, Phys. Rev. Lett. 95, 052302 (2005), arXiv:nucl-th/0504019 [nucl-th] .
Cloet et al. (2006) I. C. Cloet, W. Bentz, and A. W. Thomas, Phys. Lett. B 642, 210 (2006), arXiv:nucl-th/0605061 .
Cloët et al. (2016) I. C. Cloët, W. Bentz, and A. W. Thomas, Phys. Rev. Lett. 116, 032701 (2016), arXiv:1506.05875 [nucl-th] .
Saito et al. (2007) K. Saito, K. Tsushima, and A. W. Thomas, Prog. Part. Nucl. Phys. 58, 1 (2007), arXiv:hep-ph/0506314 .
Sammarruca (2008) F. Sammarruca, (2008), arXiv:0801.0879 [nucl-th] .
Guichon et al. (2008) P. A. M. Guichon, A. W. Thomas, and K. Tsushima, Nucl. Phys. A 814, 66 (2008), arXiv:0712.1925 [nucl-th] .
Close (1979) F. E. Close, An Introduction to Quarks and Partons (1979).
Carrillo-Serrano et al. (2016) M. E. Carrillo-Serrano, W. Bentz, I. C. Cloët, and A. W. Thomas, Phys. Lett. B759, 178 (2016), arXiv:1603.02741 [nucl-th] .
Landau (1956) L. D. Landau, Sov. Phys. JETP 3, 920 (1956).
Landau (1957) L. D. Landau, Sov. Phys. JETP 5, 101 (1957).
Landau (1959) L. D. Landau, Sov. Phys. JETP 8, 70 (1959).
Migdal (1967) A. B. Migdal, Theory of finite Fermi systems and applications to atomic nuclei (New York, Wiley, 1967).
Migdal et al. (1990) A. B. Migdal, E. E. Saperstein, M. A. Troitsky, and D. N. Voskresensky, Phys. Rept. 192, 179 (1990).
Baym and Chin (1976) G. Baym and S. A. Chin, Nucl. Phys. A262, 527 (1976).
Nozières (1964) P. Nozières, Theory of interacting Fermi systems (W.A. Benjamin, New York, 1964).
Negele and Orland (1998) J. Negele and H. Orland, Quantum Many-particle Systems (Westview Press, 1998).
Shankar (1994) R. Shankar, Rev. Mod. Phys. 66, 129 (1994), arXiv:cond-mat/9307009 .
’t Hooft (1976) G. ’t Hooft, Phys. Rev. D 14, 3432 (1976), [Erratum: Phys.Rev.D 18, 2199 (1978)].
Osipov et al. (2007) A. A. Osipov, B. Hiller, A. H. Blin, and J. da Providencia, Annals Phys. 322, 2021 (2007), arXiv:hep-ph/0607066 .
Tanimoto et al. (2020) T. Tanimoto, W. Bentz, and I. C. Cloët, Phys. Rev. C 101, 055204 (2020), arXiv:1903.06851 [nucl-th] .
Carrillo-Serrano et al. (2014) M. E. Carrillo-Serrano, I. C. Cloët, and A. W. Thomas, Phys. Rev. C 90, 064316 (2014), arXiv:1409.1653 [nucl-th] .
Saito and Thomas (1994) K. Saito and A. W. Thomas, Phys. Lett. B 327, 9 (1994), arXiv:nucl-th/9403015 .
Birse (1995) M. C. Birse, Phys. Rev. C 51, R1083 (1995), arXiv:nucl-th/9412001 .
Saito and Thomas (1995) K. Saito and A. W. Thomas, Phys. Rev. C 52, 2789 (1995), arXiv:nucl-th/9506003 .
Wallace et al. (1995) S. J. Wallace, F. Gross, and J. A. Tjon, Phys. Rev. Lett. 74, 228 (1995), arXiv:nucl-th/9407040 .
Brown and Weise (1987) G. E. Brown and W. Weise, Comments Nucl. Part. Phys. 17, 39 (1987).
Wallace (1998) S. J. Wallace, Nucl. Phys. A 631, 137C (1998).
Bentz et al. (1985) W. Bentz, A. Arima, H. Hyuga, K. Shimizu, and K. Yazaki, Nucl. Phys. A 436, 593 (1985).
Schwinger (1951) J. S. Schwinger, Phys. Rev. 82, 664 (1951).
Hellstern et al. (1997) G. Hellstern, R. Alkofer, and H. Reinhardt, Nucl. Phys. A625, 697 (1997), arXiv:hep-ph/9706551 [hep-ph] .
Gell-Mann (1962) M. Gell-Mann, Phys. Rev. 125, 1067 (1962).
Okubo (1962) S. Okubo, Prog. Theor. Phys. 27, 949 (1962).
Bentz et al. (2003) W. Bentz, T. Horikawa, N. Ishii, and A. W. Thomas, Nucl. Phys. A720, 95 (2003), arXiv:nucl-th/0210067 [nucl-th] .
Nagata and Hosaka (2004) K. Nagata and A. Hosaka, Prog. Theor. Phys. 111, 857 (2004), arXiv:hep-ph/0312161 .
Ulrych and Muther (1997) S. Ulrych and H. Muther, Phys. Rev. C 56, 1788 (1997), arXiv:nucl-th/9706030 .
Bentz and Cloët (2022) W. Bentz and I. C. Cloët, Phys. Rev. C 105, 014320 (2022), arXiv:2004.11605 [nucl-th] .
McNeil et al. (1986) J. A. McNeil, R. D. Amado, C. J. Horowitz, M. Oka, J. R. Shepard, and D. A. Sparrow, Phys. Rev. C 34, 746 (1986).
Ichii et al. (1988) S. Ichii, W. Bentz, A. Arima, and T. Suzuki, Nucl. Phys. A 487, 493 (1988).
Cohen and Furnstahl (1987) J. Cohen and R. J. Furnstahl, Phys. Rev. C 35, 2231 (1987).
Cohen (1993) J. Cohen, Phys. Rev. C 48, 1346 (1993).
Arima et al. (1987) A. Arima, K. Shimizu, W. Bentz, and H. Hyuga, Adv. Nucl. Phys. 18, 1 (1987).
Furnstahl and Serot (1987) R. J. Furnstahl and B. D. Serot, Nucl. Phys. A 468, 539 (1987).
Panda et al. (2002) P. K. Panda, M. E. Bracco, M. Chiapparini, E. Conte, and G. Krein, Phys. Rev. C 65, 065206 (2002), arXiv:nucl-th/0205051 .
Leong et al. (2023a) J. Leong, A. W. Thomas, and P. A. M. Guichon, (2023a), arXiv:2308.08987 [nucl-th] .
Okubo (1963) S. Okubo, Phys. Lett. 5, 165 (1963).
Zweig (1964) G. Zweig, “An SU(3) model for strong interaction symmetry and its breaking. Version 2,” in DEVELOPMENTS IN THE QUARK THEORY OF HADRONS. VOL. 1. 1964 - 1978, edited by D. B. Lichtenberg and S. P. Rosen (1964) pp. 22–101.
Iizuka (1966) J. Iizuka, Prog. Theor. Phys. Suppl. 37, 21 (1966).
Heiselberg and Hjorth-Jensen (2000) H. Heiselberg and M. Hjorth-Jensen, Phys. Rept. 328, 237 (2000), arXiv:nucl-th/9902033 .
Stone et al. (2021) J. R. Stone, V. Dexheimer, P. A. M. Guichon, A. W. Thomas, and S. Typel, Mon. Not. Roy. Astron. Soc. 502, 3476 (2021), arXiv:1906.11100 [nucl-th] .
Motta and Thomas (2022) T. F. Motta and A. W. Thomas, Mod. Phys. Lett. A 37, 2230001 (2022), arXiv:2201.11549 [nucl-th] .
Leong et al. (2023b) J. Leong, T. F. Motta, A. W. Thomas, and P. A. M. Guichon, Phys. Rev. C 108, 015804 (2023b).
Tolman (1939) R. C. Tolman, Phys. Rev. 55, 364 (1939).
Oppenheimer and Volkoff (1939) J. R. Oppenheimer and G. M. Volkoff, Phys. Rev. 55, 374 (1939).
Hatsuda and Kunihiro (1991) T. Hatsuda and T. Kunihiro, Z. Phys. C 51, 49 (1991).
Rehberg et al. (1996) P. Rehberg, S. P. Klevansky, and J. Hufner, Phys. Rev. C 53, 410 (1996), arXiv:hep-ph/9506436 .
Kato et al. (1993) M. Kato, W. Bentz, K. Yazaki, and K. Tanaka, Nucl. Phys. A 551, 541 (1993).
Cloët et al. (2014) I. C. Cloët, W. Bentz, and A. W. Thomas, Phys. Rev. C90, 045202 (2014), arXiv:1405.5542 [nucl-th] .