$η$-nucleus interaction from the $d+d$ reaction around the $η$ production threshold
N. Ikeno H. Nagahiro D. Jido S. Hirenzaki






                     




Department of Life and Environmental Agricultural Sciences,
Tottori University, Tottori
680-8551, Japan
Department of Physics, Nara Women's University, Nara 630-8506,
Japan
Research Center for Nuclear Physics (RCNP), Osaka University, 
Ibaraki 567-0047, Japan,
Department of Physics, Tokyo Metropolitan University, Hachioji 192-0397, Japan







The η mesic nucleus is considered to be one of the interesting exotic many body systems 
and has been studied since 1980's theoretically and experimentally.  
Recently, the formation of the η mesic nucleus in the fusion reactions
of the light nuclei such as d + d → (η + α) → X
has been proposed and the experiments have been performed by
WASA-at-COSY.
We develop a theoretical model to 
evaluate the formation rate of the η mesic nucleus in the fusion reactions and show the 
calculated results.  We find that the η bound states could be observed in the reactions in cases with the strong attractive and small absorptive η-nucleus interactions.   
We compare our results with existing data of the d + d
→η + α and the 
d + d →^3  He + N + π reactions.
We find that the analyses by our theoretical model with the
existing data can provide new information on the η-nucleus interaction.


14.40.Aqpi, K, and eta mesons36.10.GvMesonic, hyperonic and antiprotonic atoms and molecules25.60.PjFusion reactions
      
 

η-nucleus interaction from the d+d reaction around the
η production threshold
    N. Ikeno1 H. Nagahiro2,3
D. Jido4 S. Hirenzaki2 



    Received: date / Revised version: date
=============================================================================




§ INTRODUCTION



The existence of the bound states of the η meson in nucleus 
(η mesic nuclei) were predicted first by Haider and Liu <cit.>
in 1980's.
Stimulated by this theoretical result, there have been many studies of
the structure and the formation reactions of the η mesic
nucleus <cit.>.
Recently, the η-nucleon and η-nucleus interactions have
been studied theoretically in the context of the chiral symmetry of the
strong interaction and the
η mesic nucleus can be considered as one of the interesting objects
to investigate the aspects of the chiral symmetry at finite
density <cit.>.

As for the experimental studies, the first attempt to observe the η
mesic nucleus was performed by the (π^+, p) reaction with finite
momentum transfer <cit.>, 
and the interpretation of the data is still controversial <cit.>.
After that, there were many experimental searches of the bound states as
reported in
Refs. <cit.>
for example.
The systems with η in the light nucleus such as η-^3He
state also have been studied
seriously <cit.>
 and the data of the p+d →η +  ^3He reaction were studied 
to deduce η -  ^3He interaction <cit.>. 

So far, the existence of the η-Mg bound state was concluded in Ref. <cit.>.
However, we have not found any decisive evidence of the 
existence of the η bound state in lighter nuclei like He.  





Recently, the
new experiments of the d + d → (η + α) → X
reaction have been proposed and performed at WASA-at-COSY <cit.>.
In the experiments, the formation cross section of the η mesic
nucleus in the α particle in the d + d fusion reaction is
planned to be measured by observing the emitted particles 
from the decay of the η mesic nucleus 
below the η production threshold.





The formation rate of the emitting particles  
is expected to be enhanced at the resonance energy
of the η bound state formation.   
The shape of the observed spectra in Ref. <cit.>
were smooth without any clear peak
structures, and 
the upper limit of the η-nucleus formation cross section was
reported to be 
around 3–6 nb for d + d →^3He + n + π^0
reaction <cit.>.
To evaluate these upper limits, the Fermi motion of N^* in nucleus <cit.>
is also taken int account.
The upper limit for 
^3He + p + π^- 
final state is two times larger because of isospin <cit.>.








There were also the data of η production reaction d + d →η + α above the
threshold <cit.>, which
are expected to provide the valuable information on the reaction mechanism.


In this paper, we develop a theoretical model to evaluate the
formation cross section of the η bound states in the d + d
reaction and show the numerical results.   
This theoretical model can be used to deduce the
information on the η-nucleus interaction from the experimental spectra.
We explain the details of our theoretical model in section <ref>.  
We show the numerical results and compare them with the existing data
in section <ref> to deduce the information on η-nucleus
interaction, and summarize this paper in section <ref>. 





§ FORMULATION





In this section, we consider the
d + d → (η + α) → X
reaction and explain our theoretical model developed in this article.
In the experiments of this reaction at
WASA-at-COSY <cit.>, 
the total energy of the system is varied by changing the deuteron beam
momentum around the η + α threshold energy, which corresponds
to the beam momentum p=2.3 GeV/c,
and the η + α bound state production signal
is expected to be observed as peak structures of the cross section 
of the d + d →^3  He + p + π^-
and d + d →^3  He + n + π^0 reactions
in the energy region below the η-α production threshold
energy region.

Based on these considerations, we have developed a
phenomenological model described below.
In the model, the fusion and η meson production processes are
phenomenologically parameterized and the Green's function technique
is used to sum up all η-α final states.





First we formulate the transition amplitude for the d+d →η + α reaction.
We adopt the framework of the hadron reaction phenomenologically. 
We take a model in which 
the η meson production and the d+d →α fusion take place in a finite size region
as schematically pictured in Fig. <ref>.
All the information on the finiteness of
the reaction range, the spacial dimensions of the nuclear sizes, the structure
of the deuterons and the alpha and overlap of their wavefunctions is 
represented by transition form factor F(q⃗). 
We are interested in the η
production at the threshold, so that the final state, η and α, is dominated 
by s wave and, thus, the total spin-parity of the final state with a
pseudo scalar (η) and a scalar (α) bosons 
is 0^-. Deuteron having spin 1 is represented by a axial vector
boson. 
According to Lorenz invariance,
pseudoscalar 0^- state can be made 
out of two axial-vectors 1^+ by so-called anomalous coupling like
ϵ^μνρσ∂_μ A_ν∂_ρ
A_σ P S,
where A_μ and P are an axial-vector boson and a pseudoscalar
boson, respectively. S is a scalar boson.
 Thus, the interaction Hamiltonian 
may be written as

    H_ int   =    -i
       c  ϵ^ijk((∂_x_2^0∇^i_x_1 -
       ∂_x_1^0∇^i_x_2)   ϕ̂^j_d(x_1)
       ϕ̂^k_d (x_2))  
          ϕ̂_η^†(x_1) ϕ̂_α^†(x_2)   F(x_1,x_2)

where  
ϕ̂^i_d(x) is the deuteron field operator with spin index i, 
ϕ̂_η^†(x) and ϕ̂_α^†(x) are 
the creation operators for η and α, respectively, 
and c expresses the interaction strength.
The interaction strength c will be adjusted so as to reproduce the
observed cross section. 
The function F(x_1, x_2) in Eq. (<ref>) represents non-local transition
form factor of d + d →η + α, which is supposed to include the information on 
the d + d →α fusion, and 
the η meson production in the hadronic interaction such as N + N →η + N + N.
Assuming the translational invariant,

we define the Fourier transformation

    F(x_1, x_2) = ∫d^4q/(2π)^4 F(q⃗) e^i q· (x_1 - x_2).

The momentum transfer q⃗ of the reaction is large and 
all nucleons should participate in the fusion reaction equally.
Since it is hard to calculate F(q⃗) in a microscopic way,
we treat it phenomenologically and 
assume a functional form of F(q⃗) in the numerical evaluation.


Letting the wave functions of the incident deuterons labeled by d_1 and d_2 
be given by plane waves with momentum p_1 for deuteron d_1 and
p_2 for deuteron d_2
and writing the wave functions of η and α in the final state as 
ϕ_η^†(x⃗) e^-iE_η x^0 and
 ϕ_α^†(x⃗) e^-iE_α x^0
with the η and α energies E_η and E_α, respectively, 
we obtain the connected part 
of the S-matrix in the center of mass frame:

    S     =      -i
          N_d_1 N_d_2 N_η N_α 
       c ∫ d^4x_1 d^4 x_2ϵ^ijk 
      ( ∂_x_2^0∇^i_x_1 - ∂_x_1^0∇^i_x_2)
      
          [ χ_d_1^jχ_d_2^k
       e^- i p_1· x_1 e^-ip_2· x_2
       + χ_d_2^jχ_d_1^k
       e^- i p_2· x_1 e^-ip_1· x_2]
        
          ×∫d^4q/(2π)^4 F(q⃗) e^i q· (x_1 - x_2)ϕ^*_η(x⃗_1) e^iE_ηx_1^0ϕ^*_α(x⃗_2) e^iE_α x^0_2

where χ_d_1 and χ_d_2 are the spin wave functions for 
deuteron d_1 and d_2, respectively, 
and the normalization factors N_i are given as N_i = √(M_i/E_i) 
with mass M_i and energy E_i. 
Operating the derivatives onto the wave functions, we have 
the S-matrix as

    S    =     -i 2πδ(E_1 + E_2 - E_f) 
        N_d_1 N_d_2 N_η N_α
          ×
          c ∫ d x⃗_1 d x⃗_2ϵ^ijk 
       2 (E_2 p_1^i - E_1 p_2^i)χ_d_1^jχ_d_2^k
       e^ i p⃗_1·x⃗_1 e^i p⃗_2·x⃗_2
          ×∫d q⃗/(2π)^3 F(q⃗) e^-i q⃗· (x⃗_1 - x⃗_2)ϕ^*_η(x⃗_1) ϕ^*_α(x⃗_2),

where E_1 and E_2 are the energies for deuteron d_1 and d_2, 
respectively, E_f =  E_η + E_α and the integrations in terms of 
the time components provide energy conservation

    ∫ dx_1^0 dx_2^0dq_0/2π
       e^-ix_1^0(E_1-q^0 - E_η) e^-ix_2^0(E_2+q^0-E_α)
          
       = 2πδ(E_1 + E_2 -E_f).



In order to perform the spacial integrals, we introduce the relative coordinate for 
the final state, R⃗ and r⃗, defined as

    R⃗   =   m_ηx⃗_1 + M_αx⃗_2/m_η + M_α, 
    r⃗    =   x⃗_1 - x⃗_2 .

We also introduce the wave function for the relative motion of the final state, 
ϕ_f(r⃗), and assume that the center of mass motion of the η and α system is written as the plane wave. This implies that we replace the η and α wave functions as follows:

    N_η N_αϕ_η(x⃗_1) ϕ_α(x⃗_2)
      → N_f e^i p⃗_f·R⃗ϕ_f(r⃗),

with the momentum of the center of motion p⃗_f = p⃗_η + p⃗_α 
and the normalization of the wave function of relative motion ϕ_f(x⃗) 
given as

∫ dx⃗ | ϕ_f(x⃗)|^2 = 1.

In this coordinate, the S-matrix is written as

    S    =    -i  2πδ(E_1 + E_2 - E_f) 
        N_d_1 N_d_2 N_f
          ×
          c  ϵ^ijk 
       2 (E_2 p_1^i - E_1 p_2^i)χ_d_1^jχ_d_2^k
          ×∫ d R⃗ d r⃗ 
       e^ i (p⃗_1 + p⃗_2) ·R⃗  e^i (M_α/m_η+M_αp⃗_1- m_η/m_η+M_αp⃗_2)  ·r⃗
          ×∫d q⃗/(2π)^3 F(q⃗) e^-i q⃗·r⃗ϕ^*_f(r⃗) e^-i p⃗_f·R⃗  .

The integral in terms of R⃗ provides momentum conservation 

    ∫ d R⃗ e^i(p⃗_1 + p⃗_2 - p⃗_f) ·R⃗ 
       = (2π)^3δ (p⃗_1 + p⃗_2 - p⃗_f).

Introducing the Fourier transform

    ϕ_f(r⃗) = ∫d p⃗/(2π)^3ϕ̃_f (p⃗) e^-ip⃗·r⃗,

we perform the spacial integrals and obtain

    S    =    -i (2π)^4δ^(4)(p_1 + p_2 - p_f)
           N_d_1 N_d_2 N_f
          ×
          c  ϵ^ijk 
       2 (E_2 p_1^i - E_1 p_2^i)χ_d_1^jχ_d_2^k
          ×∫d q⃗/(2π)^3 F(q⃗)  
          ϕ̃^*_f (P⃗),

where P⃗ is defined as P⃗ = M_α/m_η+M_αp⃗_1-
m_η/m_η+M_αp⃗_2 - q⃗.

In the center of mass frame, E_1=E_2≡ E_d and p⃗_1 = - p⃗_2≡p⃗.  
Since the T matrix is given by S = 1 - i T (2π)^4δ^(4)(p_1 + p_2 - p_f)
 N_d_1 N_d_2 N_f, we obtain the T-matrix in the center 
of mass fram  as

    T =   4 c E_d p⃗· ( χ⃗_d_1×χ⃗_d_2) F̃(p⃗),

where we have defined F̃(p⃗) as 

    F̃(p⃗) ≡∫d q⃗/(2π)^3 F(q⃗)  
          ϕ̃^*_f( p⃗ - q⃗) .

With this T-matrix the cross section of the fusion with η production
can be obtained as

    dσ = 1/9∑_χ_d_1,χ_d_2, f|T|^2/8 p_ c.m. E_d
      (2π)^4δ^(4)(p_i - p_f) d p⃗_f/(2π)^3 2 E_f,

where p_ c.m. = | p⃗| and
we take average for the initial spin and sum up all the possible final states.
Performing the integral in terms of p⃗_f and taking spin sum, we obtain

    σ = 2π/9 c^2 p_ c.m.∑_f |F̃(p⃗)|^2δ(E_i - E_f),

where we have used E_i = 2 E_d = E_f. 


The total cross section can be written with the Green's function of the
η meson. 
Using Eq.(<ref>), we have

    ∑_f |F̃(p⃗)| ^2δ(E_f - E_i)       
       =   1/(2π)^6∑_fδ(E_f - E_i) 
       ∫ d q⃗_1  d r⃗_1 F(q⃗_1) e^-i(q⃗_1 - p⃗) ·r⃗_1ϕ^*_f(r⃗_1)
       
          ×∫ d q⃗_2  d r⃗_2 F^*(q⃗_2) e^i(q⃗_2 - p⃗) ·r⃗_2ϕ_f(r⃗_2).

The sum over the final states provides the inclusive spectrum and can be
evaluated by using Green's function method as follows:
Using the formula δ(x) = - 1/π Im1/x + iϵ 
for an infinitesimal quantity ϵ, we obtain

    ∑_fδ(E_f - E_i) ϕ_f^*(r⃗_1) ϕ_f(r⃗_2)
       =    - 1/π Im∑_fϕ^*_f(r⃗_1) 1/E_f - E_i + i ϵϕ_f(r⃗_2)    
       =    - 1/π Im∑_fϕ^*_f(r⃗_1)
        1/Ĥ - E_i + i ϵϕ_f(r⃗_2),

where we have used the fact that the wave function ϕ_f(r⃗) is
an eigenfunction of the Hamiltonian Ĥ for the η–α system.  
Equation (<ref>) is an representation of Green's operator 

(Ĥ - E_i + i ϵ)^-1 
in terms of the eigenfunction of the Hamiltonian. 
Thus, we write the Green's function in the coordinate space as 
G(E_i ; r⃗_1, r⃗_2), and we have 

    σ   =   2/9 c^2 p_ c.m.
       ×   
     (-)  Im∫ d r⃗_1 d r⃗_2
      f(r⃗_1) e^i p⃗· r_1 G(E_i; r⃗_1, r⃗_2) f^*(r⃗_2) e^-i p⃗·r⃗_2,

where we have introduced the coordinate space expression of F(q⃗) as

    f( r⃗) ≡∫dq⃗/(2π)^3 F(q⃗) e^- i q⃗·r⃗.


 
Further by assuming the spherically symmetric form of f for simplicity,
we replace f(r⃗) as,

    f(r⃗) → f(r).

And by making use of the multipole expansion of G,

    G(E_i ; r⃗_1, r⃗_2)  
    = ∑_ℓ m Y_ℓ m(r̂_1) Y_ℓ m^*(r̂_2)
    G^ℓ (E_i ; r_1, r_2),

we obtain the final form of σ as,

    σ   =    - 8π/9 c^2 p_ c.m.  Im∫ r_1^2 dr_1 r_2^2 dr_2  f(r_1) f^*(r_2)  
       ×   ∑_ℓ (2 ℓ +1)  j_ℓ(p r_1)
     G^ℓ (E_i ;  r_1,  r_2)  j_ℓ(p r_2),

where j_ℓ is the spherical Bessel function.
This expression 
is used for the numerical evaluation of the fusion cross section
of Eq. (<ref>).


We can divide the total cross section σ into two parts, 
the conversion part σ_ conv and the escape part σ_ esc 
as,

    σ = σ_ conv + σ_ esc,

based on the identity

    ImG = {G^† Im U G }+ {(1+G^†U^†)  ImG_0 (1+UG) },

where G_0 is the free Green's function of η and U the
η-nucleus potential <cit.>.
The first term of the R.H.S of Eq. (<ref>)
represents the contribution of the η meson absorption by the
imaginary part of the η-nucleus potential U and is called as the
conversion part.
The second term is the contribution from the η meson escape from
the nucleus and is called as the escape part.  

The conversion part of the cross section σ_ conv is evaluated
by the practical form written as,

    σ_ conv   =   
    - 8π/9 c^2 p_ c.m.∫ r_1^2 dr_1 r_2^2 dr_2 r_3^2 dr_3
     f(r_1) f^*(r_3) 
       ×    Im U_ opt(r_2)
      ∑_ℓ(2ℓ +1) j_ℓ(p r_1) 
       ×     
     G^ℓ * (E_i ;  r_1,  r_2)  G^ℓ (E_i ; r_2, r_3)
      j_ℓ(p r_3),

for the spherical η-α optical potential U_ opt(r).
We calculate the escape part as
σ_ esc = σ - σ_ conv.

The conversion and escape parts have the different energy dependence.
In the subthreshold energy region of the η meson production,
the total cross section σ is equal to the conversion part
σ_ conv since the energy of the η meson is insufficient
to escape from the nucleus and all η mesons must be absorbed to the
nucleus finally.
Thus, the signal of the formation of the η bound state is expected
to be observed in σ_ conv.
As shown in Eq. (<ref>), the expression of the conversion cross
section σ_ conv includes the Green's function G^ℓ
which is responsible for the peak and cusp structures in the spectrum 
as the consequences of the η-nucleus interaction such as the bound
state formation with angular momentum ℓ.
For higher partial waves ℓ without any bound states, 
the energy dependence of G^ℓ is tend to be rather mild and almost
flat in the energy region of the η production threshold.
In addition, σ_ conv also includes Im U_ opt
and, thus, the size of σ_ conv will be larger for stronger
absorptive potential.

Consequently, as shown later, the flat contribution to the
conversion spectrum gets larger in proportional to the
strength of the absorptive potential.

Above the threshold of η production,
we have also the contribution from the escape processes.
This escape part σ_ esc can be compared to the observed η
production cross section in the d + d →η + α
reaction.


We should mention here to the effects of the distortion of the initial
deuterons to the calculated results.
In Eq. (<ref>), the deuteron waves are introduced to the formula as
plane waves. The distortion effects will modify the deuteron-deuteron
relative wave function and could change the results. In the present
cases, however, we can expect the effects will be minor in the final
spectra shown in next section by the following reasons.
The energy range of the final spectra 
considered in this article is very narrow and restricted to only around
the η production threshold. Thus, the distortion effects between 
two deuterons are almost constant
in this narrow energy range and are expected to change only absolute
value of the spectra by an almost constant factor.
On the other hand, in the present analyses, we normalize the calculated
results using the experimental data of d+ d →η + α 
reaction observed above the threshold as we will see later.
Thus, the final results are expected to be insensitive to the deuteron
distortion.
We have checked qualitatively this statement
by introducing the spherical distortion factor to suppress the
contributions from the small relative coordinate region
in Eqs. (<ref>) and (<ref>),
and we have confirmed that the distortion effects to the spectra is
almost constant within the accuracy of around 15 % in the energy region
considered here.
Hence, we can neglect the deuteron distortion effects in this
article.



As for the numerical evaluation,  we assume the η-α optical potential
has the following form,

    U_ opt (r) = (V_0 + i W_0)ρ_α(r)/ρ_α(0),

where V_0 and W_0 are the parameters to determine the potential
strength at center of the α particle.
The density of the α particle
ρ_α (r) is assumed to have Gaussian form, 

    ρ_α(r) = ρ_α(0) exp[ - r^2/a^2],

with the range parameter a=1.373 fm
which reproduces the R.M.S radius of α to be 1.681 fm. 

The central density of this distribution is ρ_α(0) =
0.28  fm^-3≃ 1.6 ρ_0 with the normal nuclear density
ρ_0=0.17 fm^-3.
As the practical form of F, we assume the Gaussian as,

    F(p⃗) = (2π)^3/2(2/p_0^2 π)^3/4exp[ - p^2/p_0^2],

in this article and we treat p_0 as a phenomenological parameter.
We also show the numerical results obtained by choosing 
other functional forms for the transition form factor F(p⃗)
(f(r⃗)) in Appendix to estimate the functional form dependence
of our results.

We make a few comments on the 
difficulties for developing more microscopic model to evaluate the reaction rate. 
The first difficulty is the large momentum transfer.
We need around 1 GeV/c momentum transfer 
at the η production threshold in the center of mass frame.
The accuracy of the microscopic wave function in such high
momentum transfer region is not well investigated.
Another difficulty arises from the fact that the reaction is a fusion
reaction. In the fusion reaction, all particles in the system participate 
the reaction equally and receive large momentum transfer.
Because of these features of this reaction, we need the sufficiently 
accurate wave function of the five-body system (four nucleons and one
η meson) 
and the reliable description of the fusion and η production processes  
to perform the fully microscopic calculation.  


On the other hand, at the same time,
we can also find an advantage for the studies of
this reaction.
We can make use of the experimental data of the
d +d →η + α reaction just above the
η production
threshold <cit.>.
These data must provide us important information on the reaction 
and can be
used to fix the parameters included in the present model.
It should be also mentioned that the energy spectrum of the η-α
system is expected to be simple since the system is small
and may have only a few bound levels of η even if they exist.
The simple spectrum could be helpful to identify the bound levels from the data.






§ NUMERICAL RESULTS AND DISCUSSIONS









The theoretical model described in section <ref> includes three
parameters, which are the strength of the real and imaginary parts of
the η-α potential (V_0, W_0)
defined in Eq. (<ref>), and the parameter p_0 appeared in 
Eq.(<ref>) to determine the property of the function F which
physical meanings are explained in Eqs. (<ref>) and (<ref>).


We study first the sensitivity of the shape of the cross section to the
parameter p_0.
We show the calculated total cross section σ in Eq. (<ref>)
for p_0 = 400, 500, 600 MeV/c cases with (V_0, W_0) = -(100, 10) MeV 
as functions of the η excited energy E_η - m_η 
in Fig. <ref>.





The peak structure of the results with p_0 = 400 MeV/c is small 
and is located on the top of almost flat spectrum. 

We find that 

the structures appearing in the cross
section is insensitive to p_0 and almost same for three p_0 values.
Thus, we fix the value of this parameter to be p_0 = 500 MeV/c in the 
following numerical results
and focus on the sensitivity of the structure of the spectrum
to the η-α potential strength.




The parameter p_0 dependence of the total cross section can be
understood by considering
the change of spatial dimension of the f(r⃗) defined in
Eq. (<ref>).
For smaller p_0 value, the distribution of F in the momentum space
is more compact and that of f in coordinate space is wider.
Thus, for smaller p_0 values, we have relatively larger contributions of
higher partial wave ℓ of the η meson in
the calculation of the cross section in Eq. (<ref>), which 
are expected not to 
have any structures as function of energy around the threshold
as mentioned before.
Hence, for smaller p_0 values,
the small peak structure appears on the top of the almost flat
contribution in the spectrum.










In Fig. <ref>, we show again the calculated σ 
for the case with parameters (V_0, W_0)=-(100,10) MeV
for p_0 = 500 MeV/c to study the detail structure of the spectrum.
Three lines correspond to the total cross section σ, 
the conversion part σ_ conv,
and the escape part σ_ esc.
The η production threshold corresponds to E_η - m_η =0 and
the η-α bound states are expected to be produced in the
subthreshold region E_η - m_η < 0.
We can see in Fig. <ref> that 
the spectrum has non-trivial structures above the flat contribution whose
height seems to be about 1 in the scale of the vertical axis. 
There is a clear peak at 
E_η - m_η≃ -7 MeV which corresponds to the formation of
the η-α bound state.




The calculated escape part
is plotted with data in Fig. <ref>.
We have adjusted the height of the spectrum by the interaction strength c
in Eq. (<ref>).
The agreement of the spectrum shapes of the
calculated results and the data above the
threshold seems reasonably good in this potential parameter.
Hence, this observation implies that
this parameter set, which predicts the formation of
η-α bound state, does not contradicts to the d + d →η + α data above threshold.
We will make further comments for the comparison with the subthreshold data later
in this section.























In Fig. <ref>, 
we show the results for (V_0, W_0)=-(100,5), -(100,20) and 
-(100,40) MeV cases with p_0 = 500 MeV/c to see the effects of the
strength of the imaginary potential W_0 to the spectra.
We can see from the figures 
the structure of the spectra below threshold E_η -m_η < 0 is
sensitive to the value of the potential parameters W_0 and
is expected to be the good observables to investigate the η-nucleus
interaction.
Actually, the subthreshold spectrum with small imaginary potential W_0 =
-5 MeV clearly shows the existence of the bound state around E_η -
m_η = -7 MeV as a peak structure.

The width of the peak becomes wider and the peak height lower for larger
|W_0| value.
At the same time we can see that 
the structure of the spectra above the threshold E_η - m_η > 0
are relatively insensitive to the 
imaginary part of the η-nucleus interaction and the
value of W_0. 






We also show the calculated results for different V_0 and W_0 values
in Figs. <ref>, <ref> and
<ref>.

It could be interesting to note that the depth of the η-nucleus potential by the
chiral unitary model in Refs. <cit.>
is roughly close to -(50,40) MeV at normal nuclear density
for real and imaginary parts, respectively.
The depth of the so-called tρ potential evaluated by using the 
η-N scattering length a_η N = (0.28 + 0.19 i) fm
in Ref. <cit.> is around - (30,20) MeV at normal
nuclear density.
It should be noted that 
the parameters V_0 and W_0 adopted here in Eq. (<ref>) indicate the
potential strength at the center of the α-particle where 
ρ_α (0) ≃ 0.28 fm^-3 as defined in Eq. (<ref>).

The η-α potential also has been studied microscopically in 
Refs. <cit.>, 
where the various values of the η-α scattering length
a_ηα were reported.
We could compare our potential strength (V_0, W_0) with the scattering
length simply by the Born approximation as
a_ηα = -2 μ∫ dr r^2 U_ opt(r)
with the η-α reduced mass μ.
For example, based on this relation,
some of the potential parameters used in this
article correspond to the scattering length as,



  * (V_0, W_0)=-(100,5) MeV  ↔  a_ηα = 2.81 + 0.14 i fm

  * (V_0, W_0)=-(70,20) MeV  ↔  a_ηα = 1.97 + 0.56 i fm

  * (V_0, W_0)=-(50,40) MeV  ↔  a_ηα = 1.40 + 1.12 i fm.


It will be interesting to compare these numbers with the microscopic
scattering lengths, for example, results of the microscopic calculation
listed in Table IV in Ref. <cit.>.
We can see that the potential strengths adopted in this article are
within the range of uncertainties of microscopic calculation.




We find again the same tendencies in these results as in Fig. <ref>.
In Fig. <ref> for V_0 = -70 MeV case, we also find the bound state peak at
E_η - m_η≃ -1 MeV which is very clear for small |W_0|
case.

In the result shown in Fig. <ref> for V_0 = -50
MeV case, we find the cusp structure at the threshold energy, which
becomes less prominent for larger absorption potential.
In Fig. <ref>, 
for weaker attractive cases with V_0 = -30 MeV,
we find the step like structure at the threshold for weak absorptive
case with W_0 =-5 MeV, which becomes less prominent again for stronger
absorptive potential with larger |W_0| value.
From these figures,
we find that the total spectra around and below threshold 
E_η - m_η≤ 0 are 
sensitive to both the real and imaginary parts of the η-nucleus interaction described by
the parameters (V_0, W_0) and are good observables to obtain information on
the η-nucleus interaction.



As shown in the conversion spectra in
Figs. <ref>–<ref>,

we have the almost flat contributions in whole energy region as
mentioned in Sect. <ref>.
This flat contribution is considered to be a part of the background
cross sections of the experimental data.
Thus, we subtract the flat contribution in the conversion part in the
following numerical results to investigate the structure appearing in the
spectrum.
For this purpose, we subtract the minimum value of the conversion cross
section in the energy range shown in the following figures.




















It will be extremely interesting to show the binding energies and widths of
the η-α bound states obtained by solving the Klein-Gordon
equation using the same η-α potential used to calculate the
spectra in Figs. <ref>–<ref>.
The results are compiled in Table <ref>.

We find that only peak structures appearing in the spectra for the
strong attractive–weak absorptive potential cases correspond to the
existence of the bound states.
Other structures in the spectra may not indicate the existence of the
bound state, though they provide important information on the
η-α interaction definitely.







Now, we focus on the escape part of the spectrum which appears
above the threshold, E_η - m_η > 0, and can be compared to the experimental
data of the d + d →η + α reaction as already shown 
in Fig. <ref>
for a certain parameter set.

In Fig. <ref>, we show the calculated escape part
of the spectra for V_0 = -100 MeV cases with different
strengths of the absorptive potential together with the experimental data.
The calculated results are scaled to fit the experimental data
by changing the interaction strength c.
As we can see from Figs. <ref> and <ref>,
the shape of the escape part is relatively insensitive to the value of the potential
parameter W_0 in this case.




We also show the calculated results of the escape part with different
potential parameters in Figs. <ref>,  <ref>
 and  <ref> for different strengths of the
 attractive potential.
We find that the shape of the escape part is not very
sensitive to the V_0 and W_0 parameters and that it could be uneasy to
obtain the detail information on η-nucleus potential only from the
escape parts.
Although we have some cases which could be safely ruled out 
by these comparison, such as (V_0, W_0) = -(70,5) and -(50,5) MeV
cases, 
in which we find a distinct threshold structure 
as shown in Figs. <ref> and <ref>,
the whole shape of
the calculated results are not so much different from
that of experimental data in many cases.




















One of the best ways to obtain the decisive information on the
η-nucleus interaction may be to have the direct experimental 
observation of the η-nucleus spectrum below the threshold,
where the spectrum shape is more sensitive to the potential 
parameters as shown in Figs. <ref>–<ref>, 
especially 
if a bound state exists. 
Nevertheless, it could happen that nature might 
not 
give us 
any bound states or our experimental techniques would not
be sufficient to distinguish a less prominent peak due to
a large width. In such cases, we could deduce the information 
on the η-nucleus interaction from the absolute value of the
spectrum 
below the η-nucleus threshold 
in comparison with the η production cross section 
in the same reaction 
above the threshold. 

Here we have a model which 
can be used to 
calculate both the conversion 
and escape parts in the same footing simultaneously. 
The conversion part describes the spectrum shape induced by the
η absorption to the nucleus, 
while the escape part shows the
η production cross section. In our model, we leave 
the interaction strength of the d+d →η + α given in 
Eq. (<ref>) 
as a free parameter. We can adjust this parameter so as to 
reproduce the η production 
data 
by the escape part of our 
spectrum, and then the conversion part can be a outcome 
from the model. Since the conversion part of the spectrum
below the threshold is more sensitive to the 
η-nucleus 
interaction 
parameters, comparing the theoretical prediction and 
experimental data of the d+d reaction below the threshold,
we 
can 
deduce the information on the interaction parameters. 
























For the purpose, we show the scaled theoretical total cross section of the
d + d → (η + α) → X reaction in 
Figs. <ref>–<ref>
for various values of the parameter (V_0, W_0).
The absolute value of the cross section in these figures
are determined so that the escape part of the cross section
reproduces the experimental data of d+d →η + α 
and the nonstructural flat contributions described above
are subtracted in the conversion part. 











We should mention here that the structure of the spectra in these figures
is enhanced for smaller |W_0| values in 
Figs. <ref> and <ref>,
while it is suppressed for smaller |W_0| values in Figs. <ref>
and <ref>.
This behavior can be 
understood by considering the origin of the structure of the spectrum. 

For the strong attractive potential cases,
since the structure is dominated by the peak due to the existence of the bound state,
the peak structure becomes more prominent for the weaker imaginary
because of the smaller width of the bound state.
On the other hand, for the weak real potential cases,
since the structure is dominated by the absorptive processes,
the structure can be suppressed for the weaker absorptive potentials.



For comparison of our calculated results with experimental data
for the subthreshold energies,
we show only the conversion part, since in experiment 
the system energy is measured by observing a pion, nucleon 
and a residual nucleus emitted due to η absorption 
and this process is counted in the conversion part in the 
calculation.








In Figs. <ref>–<ref>,
we show the calculated conversion parts of the spectra which correspond to the
η absorption processes.

The obtained spectra shown in
Figs. <ref>–<ref>
can be compared to the shape of the experimental spectra on the 
background reported in 
Ref. <cit.>,
where the upper limit of the peak structure of the 
d + d →^3  He + n + π^0
 is 3–6 nb.







This implies that
the experimental upper limit for the semi-inclusive conversion 
spectrum of d+d → (η + α) → ^3He + N+
π 
including both n + π^0 and p + π^- channel
could be 
estimated to be 
3–6 nb × 3= 9–18 nb because of the isospin symmetry of the decay
channel of the η-α system.

In this case, the peak structures of the bound states in 
Figs. <ref> and <ref>
are fully rejected and strong attractive with less absorptive potentials are not
allowed.










In addition,
the upper limit provides the strong restriction to the
η-α
potential and only weak potential cases with 
small |V_0| and |W_0| values 
such as (V_0, W_0)= -(50, 5) and -(30, 5) cases
could be allowed by the limit.









In order to understand the meaning of the experimental upper limit and
the results of our analyses more clearly, we plot a contour plot of
the height of the structure appeared in the conversion spectra on the
flat contribution in the V_0 - W_0
plane in Fig. <ref>,
where the acceptable region of V_0, W_0 values can be easily
understood for each value of the upper limit of the height of the
structure in the conversion spectra
such as shown 
in Figs. <ref>–<ref>.
From this figure, the upper limit reported in
Ref. <cit.> is found to provide very valuable
information on η-nucleus interaction
and strongly suggests the small |V_0| and  |W_0| values.

However, it should be noted that the results in Fig. <ref> are
considered to be qualitative
since 
we have not considered the experimental energy resolution here.


In addition, the shapes of the structure appearing in 
Figs. <ref>–<ref>
are not simply the symmetric peak.

We also need to understand the origin of the absorptive potential and
branching ratio of various decay processes to compare our results to the
data for the specific decay mode.

Actually, we have only considered one-nucleon absorption process for
η meson here.
Multi-nucleon processes are also possible in reality <cit.>.

Thus, for deducing the quantitative information on η-nucleus interaction,
it is mandatory to make the detail comparison between 
the calculated results and data
by taking account of the realistic experimental energy resolution,
asymmetric shape of the structure appeared in the conversion spectra,
the branching ratio of the decay process of η bound states and so on,
especially for subthreshold energy region, where the spectra may have
variety of structures depending on the η-nucleus interaction
strength.



Finally, we mention the effects of the possible energy dependence of the
optical potential. The optical potential has the energy dependence in
general and the dependence could change the calculated spectra.
To simulate the energy dependence of η-α optical potential, we
adopted the energy dependence of the η-nucleon scattering length in 
Ref. <cit.> and 
we assumed the η-nucleon relative energy 
to be the quarter of that of the η-α.
The potential strength are normalized at the threshold energy.
The calculated results are shown in Figs. <ref>
and <ref>.
We have found that the energy dependence of the imaginary part of the
optical potential mainly affects the strength of the conversion part in the
spectra and changes the flat contribution of conversion part to the
slope with some gradient.
Hence, this effect could be important for the more realistic analyses.
Though, there are many theoretical models for η-nucleon scattering
length as compiled in Ref. <cit.>,
the qualitative features seems common for all models.










§ CONCLUSION

We have developed a theoretical model to evaluate the formation rate of
the η-α bound states in the d + d fusion reaction.
Because of the difficulties due to the large momentum transfer which is
unavoidable to produce η meson in the fusion reaction, we formulate
the model in a phenomenological way. 
We have shown the numerical results for the cases with the various sets
of the η-nucleus interaction parameters.

We have found that the data of η production above threshold provide
important information on the absolute strength of the reaction 
by comparing them with the escape part of the calculated results.

The upper limit of the formation cross section of the η mesic
nucleus reported in Ref. <cit.> also provides the significant
information on the strength of the η-nucleus interaction.
We would like to stress here that 
simultaneous fit to both data of d+d →η + α and 
d +d → (η + α) → X 
using our model make it possible to provide  
valuable information on η-nucleus interaction.
The results of our analyses are compiled in Fig. <ref> as
a contour plot of the V_0 - W_0 plane.










The present discussion is simply based on the value of the
upper limit of the peak structure in the fusion reaction spectrum
below the threshold. As for the further works, to make the analyses
performed in this article more quantitative and developed,
direct comparison of the spectrum shapes between the calculated
results and experiments should be necessary. For this purpose,
we should take account of experimental energy resolution in the 
calculation and consider other possibilities of the shapes of 
the spectrum structure by improving the η-nucleus optical 
potential. 


 
We acknowledge the fruitful discussion with 
P. Moskal, W. Krzemien, and M. Skurzok. 
S. H. thanks A. Gal, N. G. Kelkar, S. Wycech, E. Oset and V. Metag for fruitful comments
 and discussion in Krakow.
We also thank K. Itahashi and H. Fujioka for many discussions and
 collaborations on meson-nucleus systems.
This work was partly supported by 
JSPS KAKENHI Grant Numbers JP24540274 and JP16K05355 (S.H.), 
17K05443 (H.N.), JP15H06413 (N.I.), and JP17K05449 (D.J.)
in Japan.
 


 

 





§ APPENDIX


 In this appendix, we show the numerical results for the different
 functional form of the transition form factor.
We consider the two functions defined as

    f_1 (r⃗) =  (m/2π)^1/2e^-mr/r

with m= p_0/√(6), and 

    f_2 (r⃗) =  (λ^3/π)^1/2 e^- λ r

with λ = p_0 as different forms of the transition form factor. 
The Gaussian form defined in Eq. (<ref>) corresponds to 

    f(r⃗)   
     = (p_0^2 /2π)^3/4exp[ - p_0^2 r^2/4],

in the coordinate space.
The parameters m and λ are fixed to reproduce the same `root
mean square radius' ( ∫ |f|^2 r^2  d r⃗)^1/2 with the Gaussian form factor.

We show the calculated results with f_1(r⃗) in
Figs. <ref>–<ref>
and results with f_2(r⃗) in
Figs. <ref>–<ref>,
which correspond to Figs. <ref>,
<ref>, <ref>,  
<ref>  
obtained with f(r⃗) with p_0 =500 MeV/c.
These results are observables which can be compared
with the appropriate experimental data.
We have found that all results resemble each other
and that the numerical results are robust to the choice of the
functional form of the transition form factor.

























































































Liu



  Q. Haider and L. C. Liu,
  
  Phys. Lett. B 172, (1986) 257;
  Phys. Rev. C 34, (1986) 1845.
  





Chrien:1988gn
  R. E. Chrien et al.,
  
  Phys. Rev. Lett.  60, (1988) 2595.
  
  
  


Berger:1988ba
  J. Berger et al.,
  
  Phys. Rev. Lett.  61, (1988) 919.
  
  
  


Kohno:1989wn
  M. Kohno and H. Tanabe,
  
  Phys. Lett. B 231, (1989) 219.
  
  
  

Kohno:1990xv
  M. Kohno and H. Tanabe,
  
  Nucl. Phys. A 519, (1990) 755.
  
  
  



Chiang:1990ft
  H. C. Chiang, E. Oset and L. C. Liu,
  
  Phys. Rev. C 44, (1991) 738.
  
  
  


Sokol
G. A. Sokol and V. A. Tryasuchev, Bull. Lebedev Phys. Inst. 1991N4, 21
(1991) [Kratk. Soobshch. Fiz. 4, 23 (1991\SPLRD, 4, 21–24. 1991)]; 
G. A. Sokol, T. A. Aibergenov, A. V. Kravtsov, A. I. L'vov, and
	L. N. Pavlyuchenko, Fizika B 8, 85 (1999).


Johnson:1993zy
  J. D. Johnson et al.,
  
  Phys. Rev. C 47, (1993) 2571.
  
  
  

Waas:1997pe
  T. Waas and W. Weise,
  
  Nucl. Phys. A 625, (1997) 287.
  
  
  


Tsushima_Saito
  K. Tsushima, D. H. Lu, A. W. Thomas and K. Saito,
  
  Phys. Lett. B 443, (1998) 26;
  
  
  
  


  K. Saito, K. Tsushima, D. H. Lu and A. W. Thomas,
  
  Phys. Rev. C 59, (1999) 1203.
  
  
  
  


Hayano:1998sy
  R. S. Hayano, S. Hirenzaki and A. Gillitzer,
  
  Eur. Phys. J. A 6, (1999) 99.
  
  
  
  



Inoue:2002xw
  T. Inoue and E. Oset,
  
  Nucl. Phys. A 710, (2002) 354.
  
  
  
  

GarciaRecio:2002cu
  C. Garcia-Recio, J. Nieves, T. Inoue and E. Oset,
  
  Phys. Lett. B 550, (2002) 47.
  
  
  
  




Jido
  D. Jido, H. Nagahiro and S. Hirenzaki,
  
  Phys. Rev. C 66 (2002) 045202;
  
  
  
  

  
  
  Nucl. Phys. A 721, (2003) 665.
  
  
  
  



Nagahiro:2003iv
  H. Nagahiro, D. Jido and S. Hirenzaki,
  
  Phys. Rev. C 68, (2003) 035205.
  
  
  
  



Pfeiffer:2003zd
  M. Pfeiffer et al.,
  
  Phys. Rev. Lett.  92, (2004) 252001.
  
  
  
  

Hanhart:2004qs
  C. Hanhart,
  
  Phys. Rev. Lett.  94, (2005) 049101.
  
  
  
  



Nagahiro:2005gf
  H. Nagahiro, D. Jido and S. Hirenzaki,
  
  Nucl. Phys. A 761, (2005) 92.
  
  
  
  

Kelkar:2006zs
  N. G. Kelkar, K. P. Khemchandani and B. K. Jain,
  
  J. Phys. G 32, (2006) L19.
  
  
  
  

Jido:2008ng
  D. Jido, E. E. Kolomeitsev, H. Nagahiro and S. Hirenzaki,
  
  Nucl. Phys. A 811, (2008) 158.
  
  
  
  

Song:2008ss
  C. Y. Song, X. H. Zhong, L. Li and P. Z. Ning,
  
  Europhys. Lett.  81, (2008) 42002.
  
  
  
  




Budzanowski:2008fr
  A. Budzanowski et al. [COSY-GEM Collaboration],
  
  Phys. Rev. C 79, (2009) 012201;
  
  
  
  


  V. Jha et al. [GEM Collaboration],
  
  Int. J. Mod. Phys. A 22, (2007) 596.
  
  
  

Nagahiro:2008rj
  H. Nagahiro, D. Jido and S. Hirenzaki,
  
  Phys. Rev. C 80, (2009) 025205.
  
  
  
  


Haider:2015fea
  Q. Haider and L. C. (L. C. )Liu,
  
  Int. J. Mod. Phys. E 24 10, (2015) 1530009.
  
  
  
  





Khemchandani:2001kp
  K. P. Khemchandani, N. G. Kelkar and B. K. Jain,
  
  Nucl. Phys. A 708 (2002) 312.
  
  
  
  


Khemchandani:2003dk
  K. P. Khemchandani, N. G. Kelkar and B. K. Jain,
  
  Phys. Rev. C 68 (2003) 064610.
  
  
  
  


Khemchandani:2007ta
  K. P. Khemchandani, N. G. Kelkar and B. K. Jain,
  
  Phys. Rev. C 76 (2007) 069801.
  
  
  
  


Xie:2016zhs
  J. J. Xie, W. H. Liang, E. Oset, P. Moskal, M. Skurzok and C. Wilkin,
  
  Phys. Rev. C 95 (2017) no.1,  015202.
  
  
  
  









Krzemien:2015fsa
  W. Krzemien, P. Moskal and M. Skurzok,
  
  Acta Phys. Polon. B 46, (2015) 757.

  
  


Krzemien:2014qma
  W. Krzemien, P. Moskal and M. Skurzok,
  
  Few Body Syst.  55, (2014) 795.
  

Krzemien:2014ywa
  W. Krzemien et al. [WASA-at-COSY Collaboration],
  
  Acta Phys. Polon. B 45, (2014)  689.

  
  
  

Adlarson:2013xg
  P. Adlarson et al. [WASA-at-COSY Collaboration],
  
  Phys. Rev. C 87, (2013) 035204.

  
  


Skurzok:2016fuv
  M. Skurzok et al. [WASA-at-COSY Collaboration],
  
  Acta Phys. Polon. B 47, (2016) 503.
  
  
  
  


Skurzok:2011aa
  M. Skurzok, P. Moskal and W. Krzemien,
  
  Prog. Part. Nucl. Phys.  67, (2012) 445.
  
  
  


Adlarson:2016dme
  P. Adlarson et al.,
  
  Nucl. Phys. A 959 (2017) 102.
  
  
  
  

Kelkar:2016uwa
  N. G. Kelkar,
  
  Eur. Phys. J. A 52 (2016) no.10,  309.
  
  
  
  

Frascaria:1994va 
  R. Frascaria et al.,
  
  Phys. Rev. C 50, (1994) 537.
  
	
Willis:1997ix
  N. Willis et al.,
  
  Phys. Lett. B 406, (1997) 14.

  
	
Wronska:2005wk 
  A. Wronska et al.,
  
  Eur. Phys. J. A 26, (2005) 421.

  
	











Green

O. Morimatsu and K. Yazaki,
Nucl. Phys. A 435, (1985) 727;
A 483, (1988) 493.



	
Bhalerao:1985cr
  R. S. Bhalerao and L. C. Liu,
  
  Phys. Rev. Lett.  54 (1985) 865.
  
  
  

	








Rakityansky:1996gw
  S. A. Rakityansky, S. A. Sofianos, M. Braun, V. B. Belyaev and W. Sandhas,
  
  Phys. Rev. C 53 (1996) R2043.
  doi:10.1103/PhysRevC.53.R2043
  
  


Kelkar:2007pn
  N. G. Kelkar,
  
  Phys. Rev. Lett.  99 (2007) 210403
  doi:10.1103/PhysRevLett.99.210403
  [arXiv:0711.4066 [quant-ph]].
  
  



Fix:2017ani
  A. Fix and O. Kolesnikov,
  
  Phys. Lett. B 772 (2017) 663.
  
  
  
  


Kulpa:1998vj
  J. Kulpa and S. Wycech,
  
  Acta Phys. Polon. B 29 (1998) 3077.
  
  


Cieply:2013sya
  A. Cieplý and J. Smejkal,
  
  Nucl. Phys. A 919 (2013) 46.
  
 
  
  


Cieply:2013sga
  A. Cieplý, E. Friedman, A. Gal and J. Mareš,
  
  Nucl. Phys. A 925 (2014) 126.
  
  
  
  





