# Trapping microparticles in a structured dark focus

F. Almeida<EMAIL_ADDRESS>Departamento de Física, Pontifícia
Universidade Católica do Rio de Janeiro, 22451-900 Rio de Janeiro, RJ, Brazil
I. Sousa  Departamento de Física, Pontifícia Universidade Católica do Rio de
Janeiro, 22451-900 Rio de Janeiro, RJ, Brazil  O. Kremer  Centro de Estudos em
Telecomunicações, Pontifícia Universidade Católica do Rio de Janeiro,
22451-900 Rio de Janeiro, RJ, Brazil  B. Pinheiro da Silva  Instituto de
Física, Universidade Federal Fluminense, Niterói, Rio de Janeiro 24210-346,
Brazil  D. S. Tasca  Instituto de Física, Universidade Federal Fluminense,
Niterói, Rio de Janeiro 24210-346, Brazil  A. Z. Khoury  Instituto de Física,
Universidade Federal Fluminense, Niterói, Rio de Janeiro 24210-346, Brazil  G.
Temporão  Centro de Estudos em Telecomunicações, Pontifícia Universidade
Católica do Rio de Janeiro, 22451-900 Rio de Janeiro, RJ, Brazil  T. Guerreiro
<EMAIL_ADDRESS>Departamento de Física, Pontifícia Universidade Católica do
Rio de Janeiro, 22451-900 Rio de Janeiro, RJ, Brazil

###### Abstract

We experimentally demonstrate stable trapping and controlled manipulation of
silica microspheres in a structured optical beam consisting of a dark focus
surrounded by light in all directions - the Dark Focus Tweezer. Results from
power spectrum and potential analysis demonstrate the non-harmonicity of the
trapping potential landspace, which is reconstructed from experimental data in
agreement to Lorentz-Mie numerical simulations. Applications of the dark
tweezer in levitated optomechanics and biophysics are discussed.

††preprint: APS/123-QED

Introduction. – Light exerts forces upon matter [1]. As shown by Arthur Ashkin
[2], these forces can be used to create stable traps for nano- and microscopic
dielectric particles, with a myriad of applications from fundamental physics
[3, 4, 5, 6, 7, 8] to metrology [9, 10, 11], quantum information [12, 13] and
biology [14, 15, 16]. When the refractive index of the particle’s material is
larger than that of its surrounding medium, optical forces attract the object
towards high intensities of light. For Gaussian beam optical tweezers, the
resulting potential is approximately harmonic [17], and careful calibration of
the trap by a number of different methods [18, 19] allows for precision force
microscopy down to the molecular realm [20].

A growing interest in the fields of levitated optomechanics and optical micro-
manipulation is in enhanced force effects due to structured materials and
light beams. For example, stable Casimir trapping of refractive-index
engineered materials [21], enhanced forces in optically active nano-crystals
[22] and nitrogen-vacancy colour center ensembles [23] and probing of
structured beams using levitated nanorods [24] have been demonstrated, while
enhanced optical tweezing of meta-atoms exploiting Mie-resonances [25],
composite microspheres [26] and chiral sorting of microparticles proposed [27,
28]. Within this context of engineered nano and micro-traps, we can also
explore repulsive optical forces: in the situation that a particle has a lower
refractive index than its surrounding medium, it gets expelled from high
intensity regions of light [29]. Using structured beams [30], we can then
engineer an inverted optical trap – a dark focal region surrounded by a bright
barrier [31] – capable of trapping an object and shielding it from external
influence.

A Dark Focus Tweezer (DFT) could find many applications across physics and
biology. The optical potential generated by structured light dark traps can
have tunable non-harmonicity [32], providing a laboratory for studies of non-
linear stochastic dynamics [17] and non-Gaussian state preparation in
optomechanics [33]. Moreover, trapping objects in the dark can be extremely
beneficial in the fields of active matter and biophysics, where laser damage
limits experiments with living cells [34, 35, 36].

In this letter we report the construction of a dark optical trap for
microparticles as proposed and theoretically analysed in [32]. Earlier
experiments have employed structured light and optical bottle beams to
manipulate atoms in blue-detuned lasers [37, 38, 39, 40] and micron-sized
objects through photophoretic and thermal forces [41]. Here we demonstrate
stable trapping and controlled manipulation of microparticles in a structured
light dark focus through optical forces alone. As we will show, the DFT
induces a strongly non-harmonic potential landscape reflected as non-
Gaussianity in the statistical properties of the particle’s stochastic
trajectory. We probe the particle motion both through the its power spectrum
density as well as potential analysis and reconstruct the optical potential
landscape through matching of data with numerical simulations.

Figure 1: a) Simplified setup: a $780\text{\,}\mathrm{n}\mathrm{m}$ CW laser
is split into orthogonal polarizations by a half wave plate ($\lambda/2$) and
a polarizing beam splitter (PBS). The vertically polarized beam (Trap) is
modulated by an SLM, and directed to an objective (NA $=0.6$) to generate the
optical trap. The horizontal component (Probe) traverses the trap and is
collected by a second objective (NA $=0.4$) and used to probe the motion of a
trapped microparticle with a QPD. The trapping beam is filtered by a polarizer
(Pol). Image of the trapped particle is obtained by focusing light from a LED
onto the particle, subsequently collected by the trapping objective, filtered
by a short pass filter (SPF) and focused onto a CCD (Image). b) Behavior of
the intensity distribution around the focal point can be mimicked by varying
the relative phase $\theta$ between the Gaussian and LG modes.

The dark focus tweezer. – The DFT, sometimes also called the optical bottle
beam [31], consists of a dark focus surrounded by a bright intensity region
[32]. There are different ways of generating a DFT [42, 43, 44, 45, 46]; for
simplicity, we choose a superposition of a Gaussian ($\ell=0,p=0$) with a
Laguerre-Gauss beam ($\ell=0,p=1$) with a relative phase of $\pi$ [32]. This
choice for a bottle beam allows for an intuitive description of the optical
potential and can be readily generated using a spatial light modulator (SLM).

For a DFT of wavelength $\lambda_{0}$ in a medium of refractive index $n_{m}$,
the most important parameter is the numerical aperture $\mathrm{NA}$, from
which the beam waist $\omega_{0}=\lambda_{0}/\pi\mathrm{NA}$ and Rayleigh
range $z_{R}=n_{m}\lambda_{0}/\pi\mathrm{NA}^{2}$ can be calculated. The total
intensity of the beam is $I_{0}=2P_{0}/\pi\omega_{0}^{2}$, where $P_{0}$ is
the total beam power. Throughout this work, we adopt the intensity profile of
a DFT with $p=1$. We can also define the width $W$ and height $H$ of the
bottle as the distances between the peak values of intensity along the $x$ and
$z$ directions; these are $W=2\omega_{0},H=2z_{R}$, from which we see that the
width scales as $\mathrm{NA}^{-1}$ and height as $\mathrm{NA}^{-2}$ (see
supplemental material SM).

For a particle of radius much smaller than the beam wavelength,
$R\ll\lambda_{0}$, optical forces due to a linearly polarized light beam are
decomposed into scattering (non-conservative) and gradient (conservative)
components, both increasing with the factor
$\alpha=\left[(m^{2}-1)/(m^{2}+2)\right]$, where $m=n_{p}/n_{m}$ is the
particle-medium refractive index ratio. We are interested in situations where
$m<1$ (i.e. $n_{p}<n_{m}$) and the particles are repelled by higher
intensities of light [47]. The gradient force field can be expressed in terms
of the potential landscape,

$\displaystyle V(\vec{r})$ $\displaystyle=$ $\displaystyle-\frac{2\pi
n_{m}R^{3}}{c}\alpha\ I(\vec{r})$ (1)

where $I(\vec{r})$ is the beam intensity at position $\vec{r}$. Note that the
potential may switch from attractive to repulsive depending on the value of
$m$; although our experimental conditions do not fit the dipole regime, this
feature remains valid in our experiment. As discussed in [32], near the origin
(i.e. $\rho\ll\omega_{0}$, $z\ll z_{R}$), the potential $V(\vec{r})$ can be
expanded as a polynomial function of coordinates,

$\displaystyle V(\rho,z)\approx\frac{k_{z}}{2}z^{2}-k_{\rho
z}\rho^{2}z^{2}+\frac{k_{\rho}}{4}\rho^{4},$ (2)

where $k_{z}$ is the harmonic term strength along the axial direction and
$k_{\rho z},k_{\rho}$ denote the anharmonic potential strengths. In the
Rayleigh regime, these potential coefficients are simple functions of the beam
parameters (see SM) [32].

In our experiment trapped particles have a radius of
$R=$575\text{\,}\mathrm{n}\mathrm{m}$$, comparable to the wavelength
$\lambda_{0}=$780\text{\,}\mathrm{n}\mathrm{m}$$, a regime in which
generalized Lorentz-Mie scattering theory must be employed for the calculation
of optical forces [48]. Numerical simulation of the resulting force fields can
be performed using the toolbox presented in [49]. Eq. (2) provides a good
approximation to the potential landscape near the origin also in the
intermediate regime, with root-mean-square deviations with respect to full
numerical simulation of Lorentz-Mie scattering theory below 1% for a wide
range of particle radii (see SM).

Experimental setup. – The experimental setup for generating a dark focus
tweezer can be seen in Fig. 1a). A CW laser at
$780\text{\,}\mathrm{n}\mathrm{m}$ (Toptica DL-pro) seeds a tapered amplifier
(Toptica BoosTa) yielding $1.6\text{\,}\mathrm{W}$ of power. The beam is
divided by a half-wave plate and a polarizing beam splitter (PBS) to produce
the trapping beam and an auxiliary probe beam. The trapping beam is modulated
by an SLM (Holoeye) and sent through an objective (Olympus UPlanFLN 100x
adjustable NA $=0.6-1.3$). The resulting superposition can be monitored in a
camera providing a visualization of the beam’s transverse profile as shown in
Fig. 1b). We change the relative phase between the modes in order to mimic the
DFT’s intensity pattern along the axial direction. Image of the trapped
particle is produced by focusing light from a LED into the sample,
subsequently collecting it with the trapping objective and projecting onto a
CCD (Image camera).

SiO2 beads of radius $R=$575\text{\,}\mathrm{n}\mathrm{m}$$ (microParticles
GmbH) with refractive index $n_{p}=1.45$ are immersed in clover oil solution,
with a refractive index of $n_{m}=1.53$ and measured transmission for
$780\text{\,}\mathrm{n}\mathrm{m}$ of $\eta_{\rm clover}=85\%$. To load the
trap, we position the center of the beam at the location of a nanoparticle and
abruptly turn on the dark focus tweezer. The hydrophobic nature of the oil
increases the tendency of the silica microspheres to aggregate, occasionally
forming microdumbells in addition to single particles.

Figure 2: a)-d) Controlled SLM motion of the trapped particle. Position is
tracked by monitoring the coordinate of the brightest pixel. e)-h) Shielding
effect: a dumbbell (black contour) is repelled after approaching a trapped
particle (white contour). See supplementary videos.

Controlled particle motion and shielding. – As a first demonstration of stable
trapping in the DFT, we slightly move the beam by adjusting the SLM modulation
angle allowing for a fine control of the particle position by deflecting the
dark focus center. Fig. 2a-d) shows the iteration of four different trap
positions, with red dots marking the brightest pixel in the image
approximately corresponding to the center of the microsphere; see the
supplementary video (S1). These images are obtained by collecting the light
from an LED scattered by the particle and registered with the CCD camera, as
shown in Figure 1. By switching off the Gaussian component of the DFT
superposition and producing a pure LG mode we observed the particle is lost
from the trap.

In Gaussian optical tweezers, additional objects in the sample travelling
nearby the trapped particle are drawn into the potential landscape by the
attractive optical forces. In contrast, a particle trapped in the DFT is
shielded from the influence of these external objects due to the repulsive
optical force. Fig. 2e)-h) displays typical subsequent frames of a trapped
particle (white dotted circle) surrounded by free, passing-by particles. We
observe a microdumbbell (black dashed circle) approaching the trapped particle
and subsequently repelled by the DFT beam; see supplementary video (S2).

Power spectrum analysis. – Among the most employed techniques to calibrate
optical traps [50, 51] is the power spectrum density analysis (PSD). Analysis
of the Langevin equation for a trapped particle in a harmonic potential
reveals that the PSD has a Lorentzian form with the corner frequency parameter
proportional to the trap’s stiffness [52].

The potential associated to the DFT can be modeled by a fourth-order
polynomial in the particle’s coordinates, thus being non-harmonic (see SM).
Numerical simulations of a trapped particle in the overdamped regime subject
to quartic potentials show that the PSD of the particle motion is well fitted
by a Lorentzian function, despite the exact relation between the corner
frequency and the trap’s strength parameters being unknown beyond perturbation
theory [17]. In effect, the PSD method cannot be directly used to determine
the DFT’s strength constants, but we can use it as a consistency check between
numerical simulations of the particle motion subject to optical forces in the
intermediate regime and experimental data. This indirectly provides
information on the trap’s characteristics.

Figure 3: PSD of a particle trapped in the dark focus (red) in comparison to
background noise (grey). The effective corner frequency is
$f_{\mathrm{c,DFT}}=$ $(13.4\pm 0.7)\text{\,}\mathrm{H}\mathrm{z}$.

Due to the nature of the dark trap, scattering of photons is greatly reduced,
hindering motion detection by the traditional technique of collecting light
scattered from the trapping beam [53]. To overcome this we employ an auxiliary
weak probe beam in a Gaussian mode with polarization orthogonal to the
trapping beam. Being distinguishable and provided it has low power, the probe
beam does not significantly alter the properties of the dark trap. Moreover,
any eventual residual scattering noise due to the trapping beam can be
filtered by a suitably aligned polarizer before detection, allowing access to
the information carried by the probe alone. The probe light scattered by the
particle is collected by a second objective lens (Olympus PlanN 10x, NA =
0.25) and directed to a quadrant photodetector (QPD, New Focus 2931)
generating signals proportional to the particle’s radial and axial coordinates
(see SM for details).

Figure 4: Potential analysis. a) Position traces of the particle. b) (Top)
Kullback-Leibler divergence $D(\mathcal{P}_{\rm exp}||\mathcal{P}_{\rm sim})$
between simulation and experiment as functions of the NA. The orange region
displays the range in which simulated and experimental probability
distributions are most similar. (Bottom) Corner frequency of simulated
dynamics. The blue region displays the range where the effective corner
frequency of the simulation is within the error margin of measured value
$f_{\mathrm{c,DFT}}=$ $(13.4\pm 0.7)\text{\,}\mathrm{H}\mathrm{z}$. The
hatched area shows the intersection of both methods, NA$=0.46-0.47$. c) Fitted
normalized PDF of the centroid’s position.

We now turn to measurements of the PSD of a particle trapped in the dark
focus. The measured PSD can be seen in Fig. 9 (red dots and line) together
with the background scattering noise in the absence of a trapped particle
(grey dots and line), for comparison. A Lorentzian fit to the PSD yields an
effective corner frequency of $f_{\mathrm{c,DFT}}=$ $(13.4\pm
0.7)\text{\,}\mathrm{H}\mathrm{z}$. Numerical simulations of the trapped
particle within the DFT suggests that this value lies in the NA range between
$0.46-0.49$, where we find corner frequencies in the range between
$f_{\mathrm{c,sim}}=$ $(12.6\pm 2.3)\text{\,}\mathrm{H}\mathrm{z}$ and
$f_{\mathrm{c,sim}}=$ $(19.2\pm 7.3)\text{\,}\mathrm{H}\mathrm{z}$,
respectively. Using the scaling of width $W$ and height $H$ with the NA
introduced earlier, we estimate the order-of-magnitude of the trap size to lie
between $1.0\text{\,}\mu\mathrm{m}$ and $1.1\text{\,}\mu\mathrm{m}$.
Considering that our particles have a diameter of
$1.15\text{\,}\mu\mathrm{m}$, this suggests the view that the trapped particle
feels only loose forces within the dark focal region.

Figure 5: Reconstructed potential for centroid coordinate of a particle in the
DFT.

Potential analysis. – In thermal equilibrium and in the limit that the
conservative force dominates over dissipative forces, the position probability
density function (PDF) follows
$P(\rho,z)\propto\exp\left(-V(\rho,z)/k_{B}T\right)$ and can be reconstructed
from frames of the particle motion acquired with a CCD over long times [50].
We acquired long duration videos of a trapped particle at a rate of $15.0\
\textrm{frames}/$\text{\,}\mathrm{s}$$, from which we extract the particle’s
centroid and axial coordinates using image processing [54]. The resulting
coordinate traces can be seen in Fig. 4a). Due to the potential anharmonicity,
the position PDF is expected to be non-Gaussian; a Kolmogorov-Smirnov
hypothesis test on the position time series confirms that at $0.05$
significance.

To find the best match between the data and the quartic potential model, we
ran several simulations of the particle dynamics in the DFT parametrized by
the trap’s NA. We then extract simulated PDFs for motion along the transverse
directions and numerically compute the Kullback-Leibler (KL) divergence
between each of the simulated distributions and the marginal PDFs obtained
from the experiment. Minimizing the KL divergence between simulation and
experiment is equivalent to performing a maximum likelihood estimation of the
trap’s NA [55]. Fig. 4b) displays the KL divergence averaged over the $x$ and
$y$ directions (top plot), where we find that a potential with NA = $0.46$
best describes the measured position traces, consistent with the expectation
from Gaussian optics provided the beam waist prior to the SLM, given by NA
$\leq 0.6$. Moreover, we compute the PSD of the numerical simulations’
position data, from which we obtain the corner frequency of the Lorentzian
fit. The obtained values of corner frequency for each simulation are plotted
as a function of NA in Fig. 4b) (bottom plot), where we see that NA = $0.46$
also displays the best agreement with the PSD measurements,
$f_{\mathrm{c,sim}}=$(12.6\pm 2.3)\text{\,}\mathrm{H}\mathrm{z}$$.

From the position traces we can reconstruct the PDF, which is fitted according
to the equilibrium prediction. The resulting fit is seen in Fig. 4c). Proper
calibration of the transverse directions is achieved through independent
measurements of the CCD’s pixel size compared to a reference, while the
longitudinal direction is measured by integrating the image brightness over
the trapped sphere [53], and calibrated via comparing the root-mean-square
deviations of the matched simulation with the measured data. The reconstructed
potential at the trapped sphere’s centroid position – obtained by taking the
logarithm of the PDF – is shown in Fig. 5, with the corresponding parameters
obtained from the experiment and in comparison to full Lorentz-Mie numerical
simulation of a DFT with NA$=0.46$ in Table 1.

Parameter | Experiment | Lorentz-Mie simulation
---|---|---
$k_{z}$\text{\,}\mathrm{(}\mathrm{N}\mathrm{/}\mathrm{m}\mathrm{)}$$ | $(3.86\pm 0.06)\times 10^{-7}$ | $(2.93\pm 0.79)\times 10^{-7}$
$k_{\rho z}$\text{\,}\mathrm{(}\mathrm{N}\mathrm{/}\mathrm{{}^{3}}\mathrm{m}\mathrm{)}$$ | $(8.81\pm 0.14)\times 10^{7}$ | $(8.84\pm 0.25)\times 10^{7}$
$k_{\rho}$\text{\,}\mathrm{(}\mathrm{N}\mathrm{/}\mathrm{{}^{3}}\mathrm{m}\mathrm{)}$$ | $(2.26\pm 0.07)\times 10^{8}$ | $(1.63\pm 0.17)\times 10^{8}$
Table 1: Reconstructed potential parameters in comparison to numerical
simulations of Lorentz-Mie theory. Error bars are obtained by dividing the
experimental and simulated data into five sets and taking the standard
deviation.

Conclusion. – In summary, we have experimentally investigated a structured
light dark focus tweezer for dielectric microparticles immersed in a high
refractive index medium. We have shown stable trapping and isolation from
surrounding objects by repulsive optical forces, which induce a non-harmonic
potential landscape.

We expect the dark trap will find use both in applied and fundamental physics.
In biophysics, dark tweezers can provide stable trapping for organisms with
reduced laser heating. This is advantageous, as it has been shown that bright
tweezers hinder cell reproduction and exponentially decrease cell lifetime
even at modest trapping powers [56]. Moreover, the dark focus tweezer can be
used for vacuum optical trapping using doped nanoparticles, for instance with
rare-earth atoms [22] or Mie particles [25]. Note that implementing a vacuum
dark tweezer requires advancements in material science. Optical absorption by
particles with internal resonances typically lead to unstable dynamics and
particle loss caused by spectral and geometrical imperfections. If these
challenges can be overcome, the dark focus tweezer could provide the advantage
of a significantly reduced internal bulk temperature of the particle,
consequently reducing the decoherence effects caused by thermal emission.

###### Acknowledgements.

We acknowledge Bruno Melo, Igor Brandão, Cyril Laplane and Lukas Novotny for
helpful discussions, and Angela Duncke for help on the sample preparation
protocol. Conselho Nacional de Desenvolvimento Científico e Tecnológico
(CNPq), Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES),
Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro
(FAPERJ), Instituto Nacional de Ciência e Tecnologia de Informação Quântica
(INCT-IQ 465469/2014-0), Fundação de Amparo à Pesquisa do Estado de São Paulo
(FAPESP, processos 2021/06823-5 and 2021/06736-5), Code and data availability:
GitHub. https://github.com/QuantumAdventures/optical-bottle-beam

## References

  * Beth [1936] R. A. Beth, Phys. Rev. 50, 115 (1936).
  * Ashkin [1997] A. Ashkin, Proceedings of the National Academy of Sciences 94, 4853 (1997).
  * Delić  _et al._ [2020] U. Delić, M. Reisenbauer, K. Dare, D. Grass, V. Vuletić, N. Kiesel, and M. Aspelmeyer, Science 367, 892 (2020).
  * Windey  _et al._ [2019] D. Windey, C. Gonzalez-Ballestero, P. Maurer, L. Novotny, O. Romero-Isart, and R. Reimann, Physical review letters 122, 123601 (2019).
  * Magrini  _et al._ [2021] L. Magrini, P. Rosenzweig, C. Bach, A. Deutschmann-Olek, S. G. Hofer, S. Hong, N. Kiesel, A. Kugi, and M. Aspelmeyer, Nature 595, 373 (2021).
  * Rieser  _et al._ [2022] J. Rieser, M. A. Ciampini, H. Rudolph, N. Kiesel, K. Hornberger, B. A. Stickler, M. Aspelmeyer, and U. Delić, Science 377, 987 (2022).
  * Blakemore  _et al._ [2022] C. P. Blakemore, D. Martin, A. Fieguth, N. Priel, G. Venugopalan, A. Kawasaki, and G. Gratta, Phys. Rev. A 106, 023503 (2022).
  * Afek  _et al._ [2022] G. Afek, D. Carney, and D. C. Moore, Phys. Rev. Lett. 128, 101301 (2022).
  * Barzanjeh  _et al._ [2022] S. Barzanjeh, A. Xuereb, S. Gröblacher, M. Paternostro, C. A. Regal, and E. M. Weig, Nature Physics 18, 15 (2022).
  * Ricci  _et al._ [2022] F. Ricci, M. T. Cuairan, A. W. Schell, E. Hebestreit, R. A. Rica, N. Meyer, and R. Quidant, ACS nano 16, 8677 (2022).
  * Asano  _et al._ [2022] M. Asano, H. Yamaguchi, and H. Okamoto, Science Advances 8, eabq2502 (2022).
  * Fiaschi  _et al._ [2021] N. Fiaschi, B. Hensen, A. Wallucks, R. Benevides, J. Li, T. P. M. Alegre, and S. Gröblacher, Nature Photonics 15, 817 (2021).
  * Houhou  _et al._ [2022] O. Houhou, D. W. Moore, S. Bose, and A. Ferraro, Phys. Rev. A 105, 012610 (2022).
  * Nussenzveig [2018] H. M. Nussenzveig, European Biophysics Journal 47, 499 (2018).
  * Araújo  _et al._ [2019] G. R. d. S. Araújo, N. B. Viana, F. Gómez, B. Pontes, and S. Frases, The Cell Surface 5, 100028 (2019).
  * Bustamante  _et al._ [2021] C. J. Bustamante, Y. R. Chemla, S. Liu, and M. D. Wang, Nature Reviews Methods Primers 1, 1 (2021).
  * Suassuna  _et al._ [2021] B. Suassuna, B. Melo, and T. Guerreiro, Physical Review A 103 (2021), 10.1103/physreva.103.013110.
  * Pérez García  _et al._ [2018] L. Pérez García, J. Donlucas Pérez, G. Volpe, A. V Arzola, and G. Volpe, Nature communications 9, 1 (2018).
  * Gieseler  _et al._ [2021a] J. Gieseler, J. R. Gomez-Solano, A. Magazzù, I. P. Castillo, L. P. García, M. Gironella-Torrent, X. Viader-Godoy, F. Ritort, G. Pesce, A. V. Arzola, _et al._ , Advances in Optics and Photonics 13, 74 (2021a).
  * Bustamante and Yan [2022] C. Bustamante and S. Yan, Quarterly Reviews of Biophysics 55, e9 (2022).
  * Zhao  _et al._ [2019] R. Zhao, L. Li, S. Yang, W. Bao, Y. Xia, P. Ashby, Y. Wang, and X. Zhang, Science 364, 984 (2019).
  * Shan  _et al._ [2021] X. Shan, F. Wang, D. Wang, S. Wen, C. Chen, X. Di, P. Nie, J. Liao, Y. Liu, L. Ding, _et al._ , Nature Nanotechnology 16, 531 (2021).
  * Juan  _et al._ [2017] M. L. Juan, C. Bradac, B. Besga, M. Johnsson, G. Brennen, G. Molina-Terriza, and T. Volz, Nature Physics 13, 241 (2017).
  * Hu  _et al._ [2022] Y. Hu, J. J. Kingsley-Smith, M. Nikkhou, J. A. Sabin, F. J. Rodríguez-Fortuño, X. Xu, and J. Millen, arXiv preprint arXiv:2209.09759 (2022).
  * Lepeshov  _et al._ [2023] S. Lepeshov, N. Meyer, P. Maurer, O. Romero-Isart, and R. Quidant, Physical Review Letters 130, 233601 (2023).
  * Ali  _et al._ [2020a] R. Ali, F. A. Pinheiro, R. S. Dutra, and P. A. M. Neto, Phys. Rev. A 102, 023514 (2020a).
  * Ali  _et al._ [2020b] R. Ali, F. Pinheiro, R. Dutra, F. Rosa, and P. M. Neto, JOSA B 37, 2796 (2020b).
  * Ali  _et al._ [2021] R. Ali, R. Dutra, F. Pinheiro, and P. M. Neto, Optics Letters 46, 1640 (2021).
  * Ashkin [1970] A. Ashkin, Physical review letters 24, 156 (1970).
  * Yang  _et al._ [2021] Y. Yang, Y. Ren, M. Chen, Y. Arita, and C. Rosales-Guzmán, Advanced Photonics 3, 034001 (2021).
  * Arlt and Padgett [2000] J. Arlt and M. J. Padgett, Optics letters 25, 191 (2000).
  * Melo  _et al._ [2020a] B. Melo, I. Brandão, R. Rodrigues, A. Khoury, T. Guerreiro, _et al._ , Physical Review Applied 14, 034069 (2020a).
  * Neumeier  _et al._ [2022] L. Neumeier, M. A. Ciampini, O. Romero-Isart, M. Aspelmeyer, and N. Kiesel, arXiv preprint arXiv:2207.12539 (2022).
  * Zhang and Liu [2008] H. Zhang and K.-K. Liu, Journal of the Royal Society interface 5, 671 (2008).
  * Ehrlicher  _et al._ [2002] A. Ehrlicher, T. Betz, B. Stuhrmann, D. Koch, V. Milner, M. Raizen, and J. Käs, Proceedings of the National Academy of Sciences 99, 16024 (2002).
  * Blázquez-Castro [2019] A. Blázquez-Castro, Micromachines 10, 507 (2019).
  * Isenhower  _et al._ [2009] L. Isenhower, W. Williams, A. Dally, and M. Saffman, Optics letters 34, 1159 (2009).
  * Xu  _et al._ [2010] P. Xu, X. He, J. Wang, and M. Zhan, Optics letters 35, 2164 (2010).
  * Barredo  _et al._ [2020] D. Barredo, V. Lienhard, P. Scholl, S. de Léséleuc, T. Boulier, A. Browaeys, and T. Lahaye, Physical Review Letters 124, 023201 (2020).
  * Zhang  _et al._ [2011] P. Zhang, Z. Zhang, J. Prakash, S. Huang, D. Hernandez, M. Salazar, D. N. Christodoulides, and Z. Chen, Optics Letters 36, 1491 (2011).
  * Gong  _et al._ [2016] L. Gong, W. Liu, Q. Zhao, Y. Ren, X. Qiu, M. Zhong, and Y. Li, Scientific reports 6, 1 (2016).
  * Yelin  _et al._ [2004] D. Yelin, B. E. Bouma, and G. J. Tearney, Optics Letters 29, 661 (2004).
  * Du  _et al._ [2014] T. Du, T. Wang, and F. Wu, Optics Communications 317, 24 (2014).
  * Wei  _et al._ [2005] M.-D. Wei, W.-L. Shiao, and Y.-T. Lin, Optics Communications 248, 7 (2005).
  * Yang  _et al._ [2019] Z. Yang, X. Lin, H. Zhang, X. Ma, Y. Zou, L. Xu, Y. Xu, and L. Jin, Applied Optics 58, 2471 (2019).
  * Whyte and Courtial [2005] G. Whyte and J. Courtial, New Journal of Physics 7, 117 (2005).
  * Grynberg  _et al._ [2010] G. Grynberg, A. Aspect, and C. Fabre,  _Introduction to quantum optics: from the semi-classical approach to quantized light_ (Cambridge university press, 2010).
  * Jones  _et al._ [2015] P. Jones, O. Maragó, and G. Volpe,  _Optical tweezers_ (Cambridge University Press Cambridge, UK:, 2015).
  * Nieminen  _et al._ [2007] T. A. Nieminen, V. L. Loke, A. B. Stilgoe, G. Knöner, A. M. Brańczyk, N. R. Heckenberg, and H. Rubinsztein-Dunlop, Journal of Optics A: Pure and Applied Optics 9, S196 (2007).
  * Gieseler  _et al._ [2021b] J. Gieseler, J. R. Gomez-Solano, A. Magazzù, I. P. Castillo, L. P. García, M. Gironella-Torrent, X. Viader-Godoy, F. Ritort, G. Pesce, A. V. Arzola, K. Volke-Sepúlveda, and G. Volpe, Advances in Optics and Photonics 13, 74 (2021b).
  * Melo  _et al._ [2020b] B. Melo, F. Almeida, G. Temporão, and T. Guerreiro, Optics Express 28, 16256 (2020b).
  * Berg-Sørensen and Flyvbjerg [2004] K. Berg-Sørensen and H. Flyvbjerg, Review of Scientific Instruments 75, 594 (2004).
  * Tebbenjohanns  _et al._ [2019] F. Tebbenjohanns, M. Frimmer, and L. Novotny, Physical Review A 100, 043821 (2019).
  * Forsyth and Ponce [2002] D. A. Forsyth and J. Ponce,  _Computer vision: a modern approach_ (prentice hall professional technical reference, 2002).
  * Rice [2003] J. A. Rice,  _Mathematical statistics and data analysis_ (China machine press, 2003).
  * Pilát  _et al._ [2017] Z. Pilát, A. Jonáš, J. Ježek, and P. Zemánek, Sensors 17, 2640 (2017).

## Appendix A: The Dark Focus Tweezer

The Dark Focus Tweezer (also known as optical bottle beam) used in this work
consists in a superposition of a Gaussian and a Laguerre-Gauss (LG) mode with
$\ell=0,p\neq 0$ with a relative phase of $\pi$. The LG mode functions read,

$\displaystyle
u_{\ell,p}(\rho,\phi,z)=\sqrt{\frac{2}{\pi\omega(z)^{2}}}\sqrt{\frac{p!}{(|\ell|+p)!}}\times$
(3)
$\displaystyle\left(\frac{\sqrt{2}\rho}{\omega(z)}\right)^{|\ell|}\mathrm{L^{|\ell|}_{p}}\left(\frac{2\rho^{2}}{\omega(z)^{2}}\right)\exp\left[-\frac{\rho^{2}}{\omega(z)^{2}}\right]\times\exp[ik_{m}z+ik_{m}\frac{\rho^{2}}{2R(z)}-i\zeta(z)+i\ell\phi],$

where $k_{m}$ is the wavenumber in the medium and $\omega(z)$, $R(z)$,
$\zeta(z)$ and $\mathrm{L^{|\ell|}_{p}}$ are the beam width, the wavefront
radius, the Gouy phase and the Associated Laguerre polynomial. These
quantities are,

$\displaystyle\omega(z)$ $\displaystyle=$
$\displaystyle\omega_{0}\sqrt{1+\frac{z^{2}}{z^{2}_{R}}};$ (4) $\displaystyle
R(z)$ $\displaystyle=$ $\displaystyle
z\left(1+\frac{z_{R}^{2}}{z^{2}}\right);$ (5) $\displaystyle\zeta(z)$
$\displaystyle=$ $\displaystyle(2p+|\ell|+1)\arctan\frac{z}{z_{R}};$ (6)
$\displaystyle\mathrm{L_{p}^{|\ell|}}(x)$ $\displaystyle=$
$\displaystyle\sum^{p}_{i=0}\frac{1}{i!}\binom{p+|\ell|}{p-i}(-x)^{i}$ (7)

where the Rayleigh range ($z_{R}$) and the beam waist ($\omega_{0}$) are
defined as

$\displaystyle\omega_{0}=\frac{\lambda_{0}}{\pi\mathrm{NA}}\,,\quad
z_{R}=\frac{n_{m}\lambda_{0}}{\pi\mathrm{NA}^{2}}$ (8)

with $\lambda_{0}$ being the wavelength in vacuum, $n_{m}$ the medium
refractive index and $\mathrm{NA}$ the numerical aperture. Throughout this
work we consider linearly polarized fields.

The intensity profile of a Gaussian beam with total power $P_{0}$ is simply
$I_{G}(\mathbf{r})=P_{0}\,|u_{0,0}(\mathbf{r})|^{2}$. The resulting intensity
profile of the DFT beam with the same total power $P_{0}$ reads,

$\displaystyle I_{\mathrm{B,p}}(\rho,z)$ $\displaystyle=$ $\displaystyle
P_{0}\,\left|\frac{u_{0,0}(\mathbf{r})-u_{0,1}(\mathbf{r})}{\sqrt{2}}\right|^{2}$
$\displaystyle=$
$\displaystyle\frac{P_{0}}{\pi\omega(z)^{2}}\exp\left[-\frac{2\rho^{2}}{\omega(z)^{2}}\right]\times\bigg{[}1-2\cos\left(2p\arctan\frac{z}{z_{R}}\right)\mathrm{L^{0}_{p}}\left(\frac{2\rho^{2}}{\omega(z)^{2}}\right)+\
\mathrm{L^{0}_{p}}\left(\frac{2\rho^{2}}{\omega(z)^{2}}\right)^{2}\bigg{]}\,,$

where $P_{0}$ is the total power of the beam. Fig. 6 shows the resulting
normalized intensity profile of a DFT with $p=1$, which we adopt throughout
this work.

Figure 6: Intensity landscape of the dark focus tweezer generated by
superposing a Gaussian and Laguerre-Gauss mode with $\ell=0,p=1$ with a
relative phase of $\pi$. The beam intensity is normalized by
$I_{0}=P_{0}/(\pi\omega_{0}^{2})\,$.

Following [51], we define the width $W$ and height $H$ of the OBB as the
distance between the two intensity maxima surrounding the dark focus along the
$x$ axis ($z$ axis). These can be found by solving

$\displaystyle dI_{\mathrm{B,p}}(x,0,0)/dx|_{x=W/2}=0\,,\ \ \ \ \
dI_{\mathrm{B,p}}(z,0,0)/dz|_{z=H/2}=0\,,$ (10)

yielding $W=2\omega_{0},H=2z_{R}$ for $p=1$. Note the width of the OBB scales
as $\mathrm{NA}^{-1}$, while the height scales as $\mathrm{NA}^{-2}$.

The DFT intensity in the Rayleigh regime ($R\ll\lambda_{0}$) gives rise to a
scattering (non-conservative) and gradient (conservative) force fields,
respectively given by [48],

$\displaystyle\vec{F}_{{\mathrm{B,p}}}^{(s)}(\vec{r})$ $\displaystyle=$
$\displaystyle\hat{z}\frac{128\pi^{5}R^{6}}{3c\lambda_{0}^{4}}\left(\frac{m^{2}-1}{m^{2}+2}\right)^{2}n_{m}^{5}I_{p}(\vec{r})$
(11) $\displaystyle\vec{F}_{{\mathrm{B,p}}}^{(g)}(\vec{r})$ $\displaystyle=$
$\displaystyle\frac{2\pi
n_{m}R^{3}}{c}\left(\frac{m^{2}-1}{m^{2}+2}\right)\nabla I_{p}(\vec{r})$ (12)

where $m=n_{p}/n_{m}$ is the particle-medium refractive index ratio and
$I_{{\mathrm{B,p}}}(\vec{r})$ is the intensity distribution of the beam.

The potential landscape associated to the gradient force near the origin is
given by an approximate polynomial potential

$\displaystyle\frac{V_{{\mathrm{B,p}}}(\rho,z)}{V_{0}}\approx\mu^{\prime}z^{2}-\eta^{\prime}\rho^{2}z^{2}+\chi^{\prime}\rho^{4}$
(13)

with,

$\displaystyle\mu^{\prime}=\frac{4p^{2}}{z_{R}^{2}}\ ,\
\eta^{\prime}=\frac{8p^{2}(p+1)}{\omega_{0}^{2}z_{R}^{2}}\ ,\
\chi^{\prime}=\frac{4p^{2}}{\omega_{0}^{4}}$ (14)

and $V_{0}=[2\pi n_{m}R^{3}(m^{2}-1)/c(m^{2}+2)]I_{0}$. Numerical simulations
of Lorentz-Mie scattering theory validate that in the intermediate regime, the
same potential form (13) can be considered, despite the coefficients
$\mu^{\prime},\eta^{\prime},\chi^{\prime}$ are no longer expressed as simple
functions of the beam parameters.

As explained in the main text, we generate the OBB using an SLM. Figure 7
depicts the beam profile at the focal plane, together with a corresponding fit
of the beam intensity according to (Appendix A: The Dark Focus Tweezer). We
observe excellent agreement between the expected beam profile and the mode
generated by the SLM.

Figure 7: Intensity plot along the $x$ direction of OBB transverse profile at
the dark focus ($\theta=\pi$), with corresponding fit. Considering the optical
setup (telescope after the SLM), a profile beam waist of
$w_{\mathrm{expected}}=$250\text{\,}\mu\mathrm{m}$$ is expected. The beam
radius obtained by fitting (Appendix A: The Dark Focus Tweezer) to the
transverse profile along the $x$ direction is $w_{\rm
fit}=$251(3)\text{\,}\mu\mathrm{m}$$, demonstrating the generation of a good
quality OBB.

## Appendix B: Effective Potential Validation

As discussed in the main text, the potential in the dipole regime can be
expanded as a polynomial function of $\rho$ and $z$. From Eq. (13), it is
possible to find expressions for the forces in the $x$, $y$ and $z$ direction.
To validate the quartic polynomial model in the intermediate regime, we
simulated the optical forces considering particles of different radii using
the MatLab toolbox described in [49]. We fit a polynomial to the force vectors
$f_{x}$, $f_{y}$ and $f_{z}$ obtained from the simulation, returning
estimations for each of the components, namely $\hat{f}_{x}$, $\hat{f}_{y}$
and $\hat{f}_{z}$. The quality of the fit in each axis can be evaluated using
a root-mean-squared error (RMSE) divided by the root-mean-square force.
Finally an average RMSE over each axis is considered,

$\mathrm{RMSE_{avg}}=\frac{1}{3}\sum_{i\in\\{x,y,z\\}}\Bigg{\\{}\sqrt{\frac{\sum\limits_{j=1}^{N}(f_{i}^{j}-\hat{f}_{i}^{j})^{2}}{\sum\limits_{j=1}^{N}f_{i}^{j\,2}}}\Bigg{\\}},$
(15)

where $N$ represents the number of points used for discretization during
simulation.

Fig. 8 shows the $\mathrm{RMSE_{avg}}$ as a function of particle radius $R$.
Different values of NA were considered to ensure the approximation is valid
under variations of the trapping beam focusing. The maximum error encountered
is always less than $1.50\%$, for a particle radius of
$350\text{\,}\mathrm{n}\mathrm{m}$ and an NA of $0.58$. For the experimental
conditions described in the main text, the error is $0.416\%$, validating the
quartic potential model within the experimental parameters.

Figure 8: $\mathrm{RMSE_{avg}}$ of the polynomial fit for the optical forces
as a function of the particle radius. Simulations were executed considering NA
= 0.40 ( ), NA = 0.46 ( ), NA = 0.52 ( ) and NA = 0.58 ( ). The conditions in
which the experiment was conducted are represented by the black circle
($\bullet$).

## Appendix C: Probe beam calibration

Figure 9: Probe beam calibration: PSD of probe beam scattering at different
intensities for a particle trapped at fixed power $P_{T}=$
$39\text{\,}\mathrm{m}\mathrm{W}$, with corner frequency $f_{c}=$ $(428.8\pm
14)\text{\,}\mathrm{H}\mathrm{z}$ (red). The probe beam power $P_{P}^{i}$
($i=1,2,3$) is progressively decreased: $P_{P}^{1}=$
$91\text{\,}\mathrm{m}\mathrm{W}$, $f_{c}=$ $(1272.2\pm
40)\text{\,}\mathrm{H}\mathrm{z}$ (blue); $P_{P}^{2}=$
$45\text{\,}\mathrm{m}\mathrm{W}$, $f_{c}=$$(770.8\pm
19)\text{\,}\mathrm{H}\mathrm{z}$ (orange); $P_{P}^{3}$ =
$19\text{\,}\mathrm{m}\mathrm{W}$, $f_{c}=$ $(492.0\pm
22)\text{\,}\mathrm{H}\mathrm{z}$ (green). For comparison, the PSD of the
probe beam in the absence of a trap is also shown (grey). Probe beam powers
below $19\text{\,}\mathrm{m}\mathrm{W}$ allows for position read-out without
significant disturbance to the trap. Each PSD curve is normalized by its
corresponding plateau value $\mathrm{P_{xx}^{max,i}}$ obtained from the
Lorentzian fit.

To determine the proper alignment and optimal power of the probe beam we
perform standard Gaussian tweezer experiments at various probe beam powers
with a particle immersed in aqueous solution at a fixed trap power of
$P_{T}=$39\text{\,}\mathrm{m}\mathrm{W}$$. As the probe power is decreased,
the measured PSD approaches that of the trapping beam alone. The power at
which the probe and trap PSDs are indistinguishable determines the power at
which the probe beam causes negligible influence on the trap within
experimental uncertainties; this power was measured to be
$P_{P}\leq$19\text{\,}\mathrm{m}\mathrm{W}$$. The result of this calibration
measurement is shown in Fig. 9a). See the caption for details.

## Appendix D: Absorbed power in the dark focus and trap performance

To quantify the energy absorption in the dark focus in comparison to a
standard bright tweezer we turn to the formula for the absorbed power in terms
of the beam’s and particle’s parameters [10],

$\displaystyle P_{\mathrm{abs}}=12\pi
I_{0}\frac{V}{\lambda}\mathfrak{I}\left(\frac{\epsilon-1}{\epsilon+2}\right)$
(16)

where $I_{0}$ is the light intensity and $A_{p}\equiv V/\lambda$ is the
effective cross-section of the particle, so that $I_{0}V/\lambda$ represents
the optical power effectively seen by the particle. A meaningful comparison
between the usual Gaussian bright tweezer and the proposed dark focus tweezer
can be built from the following absorption ratio,

$\displaystyle\eta_{\mathrm{abs}}=\frac{P^{\,\prime}_{\mathrm{abs}}}{P_{\mathrm{abs}}}=\frac{P_{B}}{P_{G}}\times\frac{\int_{A_{p}}|u_{B}(\mathbf{r})|^{2}d^{2}\mathbf{r}}{\int_{A_{p}}|u_{G}(\mathbf{r})|^{2}d^{2}\mathbf{r}}$
(17)

where $u_{G}$ is the Gaussian field distribution used in a standard bright
tweezer, $u_{B}$ is the optical bottle beam (OBB) field distribution used in
the dark focus tweezer and $P_{G},P_{B}$ are the powers in the Gaussian and
dark focus traps, respectively. Note that the only free parameters in the
definitions of the modes $u_{B}$ and $u_{G}$ are the total power and the beam
waist $w_{0}$, and both $u_{B}$ and $u_{G}$ are normalized with respect to
integration over the whole transverse plane, while their intensity
distributions are completely different.

In order to give a quantitative estimation of the ratio
$\eta_{\mathrm{abs}}\,$, we assume that the absolute value of the
polarizability of the particle immersed in the surrounding medium is the same
in both the dark focus and Gaussian traps. The radius of the Silica particles
used in our experiment is $R=$575\text{\,}\mathrm{n}\mathrm{m}$$ and the laser
wavelength is $\lambda=$780\text{\,}\mathrm{n}\mathrm{m}$\,$, giving
$A_{p}=1.02\times 10^{-12}\ \mathrm{m}^{2}$.

If we set the waists and total power of the dark focus and Gaussian tweezers
to be the same, that is $P_{B}=P_{G}$, using the above parameters for particle
radius, wavelength $\lambda$, cross-section $A_{P}$ and the mode functions
employed in the experiment, we find an absorption ratio of
$\eta_{\mathrm{abs}}=0.045$. Moreover, the value of $\eta_{\mathrm{abs}}$
decreases with the particle radius, as expected for the dark focus tweezer.
The numerical value of $\eta_{\mathrm{abs}}$ as a function of effective
particle radius $R_{\mathrm{eff}}=\sqrt{A_{p}/4\pi}$ is indicated in Figure
10.

Figure 10: Ratio $\eta_{\mathrm{abs}}$ as a function of the effective particle
radius $R_{\mathrm{eff}}$.

It is interesting to compare the dark focus trapping characteristics to that
of a standard Gaussian tweezer. Note that the trapping potential associated
with $u_{B}$ has a quartic dependence with the radial coordinate
($V\sim\frac{k_{\rho}}{4}\rho^{4}$) over the transverse plane, so it is not
suitable to compare $k_{\rho}$ with the spring constant of the quadratic
potential produced by $u_{G}$. Once the powers $P_{G}$, $P_{B}$ and beam waist
$w_{0}$ are fixed for both beams, however, a comparison between the trapping
potentials can be made through the corresponding trap depths.

The reduction in absorption in the dark focus trap compared to the absorption
in a Gaussian tweezer of the same power and waist is accompanied by a
reduction in trap depth along the transverse plane, as can be seen in the Left
plot in Figure 11. A straightforward calculation using the mode functions
shows that,

$\displaystyle V_{0,B}=\frac{2}{e^{2}}V_{0,G}\approx 0.27V_{0,G}$ (18)

where $V_{0,B},V_{0,G}$ are the trap depths in the dark focus and Gaussian
tweezers, respectively. Therefore, we see that a significant 20-fold reduction
in absorbed power ($\eta_{\mathrm{abs}}=0.045$) comes at the expense of a
4-fold reduction in trap depth. Note that in our experiment the observed trap
depth was on the order of $V_{0,B}\approx 100k_{B}T$, guaranteeing that the
particle was well within the trapping potential.

Figure 11: Comparison between the dark focus and Gaussian tweezers for the
transverse directions and corresponding absorption ratios. Left: For the same
optical power and waist, the dark focus presents an absorption ratio of
$\eta_{\mathrm{abs}}=0.045$, i.e. the particle absorbs about $4.5\%$ of the
power absorbed in a standard optical trap. This reduction in absorption comes
at the expense of a reduction in trap depth of $V_{0,B}\approx 0.27V_{0,G}$.
Right: Increasing the dark focus power to about a factor of 4 the trap depth
of both tweezers can be matched, yielding equivalent trap performances with an
absorption ratio of $\eta_{\mathrm{abs}}\approx 0.17$. This is a substantial
reduction in absorption, i.e. $\approx 17\%$ of the total absorption in a
standard Gaussian trap. In both plots, the zero potential energy reference was
set so that the two minima coincide.

Alternatively, we can set the modes $u_{B}$ and $u_{G}$ at different powers
such that their trap depths are equal, as shown in the Right plot in Figure
11. In this case, the power in the dark focus has to be approximately 4 times
larger than in the Gaussian beam, resulting in an approximately 4-fold
increase in absorption. This yields an absorption ratio of
$\eta_{\mathrm{abs}}\approx 0.17$, showing that the trap performance of a
standard Gaussian beam can be matched by the dark focus tweezer while at the
same time achieving a substantial reduction in absorbed power, i.e. $\approx
17\%$ of the total absorption in a standard Gaussian trap.

Figure 12: Comparison between the longitudinal trapping potentials in the dark
focus ($V_{\mathrm{B}}$) and Gaussian tweezer ($V_{\mathrm{G}}$) with equal
waists. Left: The two types of trapping beams are assumed to have equal power.
The dark focus presents a reduction of about $50\%$ in trap depth, however,
near the origin it provides a quadratic potential twice stiffer than a
Gaussian trap. Right: The dark focus beam power is assumed to be about 4 times
larger than the Gaussian beam, in order to match the transverse potential
depths. Besides being stiffer, the dark focus longitudinal trapping potential
is approximately twice deeper than the Gaussian one.

Concerning the longitudinal direction, the trapping potentials of both the
bottle and Gaussian beams have a quadratic dependence with the $z$ coordinate
to leading order, so the spring constants can be directly compared. As shown
in the Left plot in Figure 12, a dark focus with the same power as the
Gaussian tweezer presents a reduction of about $50\%$ in trap depth, however,
presenting a stiffer potential near the origin. From Eq. (13) of [32], the
dark focus trap longitudinal spring constant is twice larger than that of a
Gaussian beam of equal power. Increasing the power of the dark focus by a
factor of four to match the transverse trapping characteristics of a Gaussian
beam, as discussed above and shown in the Right plot in Figure 12, the
longitudinal potential trap depth in the dark focus becomes approximately
twice that of the Gaussian.


# Euler sums of generalized harmonic numbers and connected extensions

Mümün CAN, Levent KARGIN, Ayhan DİL, Gültekin SOYLU
Department of Mathematics, Akdeniz University, Antalya, Turkey

###### Abstract

This paper presents the evaluation of the Euler sums of generalized
hyperharmonic numbers $H_{n}^{\left(p,q\right)}$

$\zeta_{H^{\left(p,q\right)}}\left(r\right)=\sum\limits_{n=1}^{\infty}\dfrac{H_{n}^{\left(p,q\right)}}{n^{r}}$

in terms of the famous Euler sums of generalized harmonic numbers. Moreover,
several infinite series, whose terms consist of certain harmonic numbers and
reciprocal binomial coefficients, are evaluated in terms of Riemann zeta
values.

##  1 Introduction

The classical Euler sum $\zeta_{H}\left(r\right)$ is the following Dirichlet
series

$\zeta_{H}\left(r\right)=\sum\limits_{n=1}^{\infty}\dfrac{H_{n}}{n^{r}},$

where $H_{n}$ is the $n$th harmonic number. This series is also known as the
harmonic zeta function. The famous Euler’s identity for this sum is [14, 22,
30]

$2\zeta_{H}\left(r\right)=\left(r+2\right)\zeta\left(r+1\right)-\sum_{j=1}^{r-2}\zeta\left(r-j\right)\zeta\left(j+1\right),\text{
}r\in\mathbb{N}\backslash\left\\{1\right\\},$ (1)

where $\zeta\left(r\right)$ is the classical Riemann zeta function (for more
details, see for instance [40]). Many generalizations of Euler sums (the so
called Euler-type sums) are given using generalizations of harmonic numbers
(see [2, 4, 3, 7, 8, 10, 21, 26, 37, 36, 41, 43, 45, 46, 47, 48, 49]).
Evaluation of Euler-type sums and construction of closed forms are active
fields of study in analytical number theory. Furthermore [3, 7, 10, 16, 17]
are some of the studies that make this area interesting in the sense that
Euler sums have potential applications in quantum field theory and knot
theory, especially in evaluation of Feynman diagrams.

Euler actually considered also the more general form [5, 19, 22, 30]

$\zeta_{H^{\left(p\right)}}\left(m\right)=\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{n^{m}},$
(2)

where $H_{n}^{\left(p\right)}$ defined by

$H_{n}^{\left(p\right)}=1+\frac{1}{2^{p}}+\frac{1}{3^{p}}+\cdots+\frac{1}{n^{p}},\text{
}\left(p\in\mathbb{Z},n\in\mathbb{N}\right),$

is the $n$th partial sum of $\zeta\left(p\right)$ and is called the $n$th
generalized harmonic number for $p>1$. In particular,
$H_{0}^{\left(p\right)}=0$ and $H_{n}^{\left(1\right)}=H_{n}$, the $n$th
harmonic number. When $p\leq 0$ it is called sum of powers of integers.

One of the most important issues here is to write Euler-type sums as
combinations of the Riemann zeta function as in (1). This problem has remained
important for various Euler-type sums from the era of Euler to the present
day. It’s shown by Euler himself that, the cases of $p=1$, $p=q$, $p+q$ odd,
and for special pairs $(p,q)\in\\{(2,4),(4,2)\\}$, the sums of the form (2)
have evaluations in terms of the Riemann zeta function (see [5, 19, 22, 30]).
There is a very comprehensive literature on this subject, both theoretical and
numerical ([1, 6, 7, 10, 19, 14, 20, 21, 29, 36, 41, 42, 43, 44, 45, 48]). One
of these results; the Euler identity (1) was further extended in the works of
Borwein et al. [6] and Huard et al. [24]. For odd weight $N\geq 3$ and
$p=1,2,\ldots,N-2,$ we have [24, Theorem 1] (or [6, p. 278])

$\displaystyle\zeta_{H^{\left(p\right)}}\left(N-p\right)$
$\displaystyle=\left(-1\right)^{p}\sum_{j=0}^{\left[\left(N-p-1\right)/2\right]}\binom{N-2j-1}{p-1}\zeta\left(N-2j\right)\zeta\left(2j\right)$
(3)
$\displaystyle+\left(-1\right)^{p}\sum_{j=0}^{\left[p/2\right]}\binom{N-2j-1}{N-p-1}\zeta\left(N-2j\right)\zeta\left(2j\right)-\zeta\left(0\right)\zeta\left(N\right).$

Moreover, these so called ”linear Euler sums” satisfy a simple reflection
formula

$\zeta_{H^{\left(p\right)}}\left(r\right)+\zeta_{H^{\left(r\right)}}\left(p\right)=\zeta\left(p+r\right)+\zeta\left(p\right)\zeta\left(r\right).$
(4)

Considering nested partial sums of the harmonic numbers, Conway and Guy [18]
introduced hyperharmonic numbers for an integer $r>1$ as

$h_{n}^{(r)}=\sum_{k=1}^{n}h_{k}^{(r-1)},\text{ }n\in\mathbb{N},$

with $h_{n}^{(0)}=1/n,$ $h_{n}^{(1)}=H_{n}$ and $h_{0}^{(r)}=0.$ Hyperharmonic
numbers are also important because they build a step in the transition to the
multiple zeta functions (see [26, 42]). Dil and Boyadzhiev [20] extended the
Euler’s identity (1) to the Euler sums of the hyperharmonic numbers:

$\zeta_{h^{\left(q\right)}}\left(r\right)=\sum_{n=1}^{\infty}\frac{h_{n}^{\left(q\right)}}{n^{r}},\text{
}(r>q),$ (5)

as

$\displaystyle\zeta_{h^{\left(q\right)}}\left(r\right)$
$\displaystyle=\frac{1}{\left(q-1\right)!}\sum_{k=1}^{q}\genfrac{[}{]}{0.0pt}{}{q}{k}$
(6)
$\displaystyle\times\left\\{\zeta_{H}\left(r-k+1\right)-H_{q-1}\zeta\left(r-k+1\right)+\sum_{j=1}^{q-1}\mu\left(r-k+1,j\right)\right\\},$

where $\genfrac{[}{]}{0.0pt}{}{q}{k}$ is the Stirling number of the first kind
and

$\mu\left(r,j\right)=\sum_{n=1}^{\infty}\frac{1}{n^{r}\left(n+j\right)}=\sum_{k=1}^{r-1}\frac{\left(-1\right)^{k-1}}{j^{k}}\zeta\left(r+1-k\right)+\left(-1\right)^{r-1}\frac{H_{j}}{j^{r}}.$
(7)

Formula (6) was the general form of the results obtained for some special
values of $q$ and $r$ in the study of [29].

Studies on evaluating Euler sums (2) and (5) in terms of Riemann zeta values
$\zeta\left(m\right)$ have motivated researchers to find representations
harmonic number series of the forms

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{\left(n+m\right)^{r}\dbinom{n+m+l}{l}},\text{
}\sum_{n=1}^{\infty}\frac{h_{n}^{\left(q\right)}}{n\dbinom{n+q}{q}}.$

It has been shown that some families of these type of series can be evaluated
in terms of Euler sums and Riemann zeta values (see for example for $m=0,$
$p=1,$ $r\in\left\\{0,1\right\\}$ [32, 34, 38], for $m=r=0$ [33], for $m=0$
[32, 44], for $m>0,$ $r=1,$ $p\in\left\\{1,2\right\\}$ [35, 39] and for the
series involving hyperharmonic numbers [20]). We would like to emphasize that
in some studies these type of series have been expressed in terms of
hypergeometric series [13, 14, 29, 35, 36, 39].

In this work we mainly concentrate on generalized hyperharmonic numbers
defined as (see [21])

$H_{n}^{\left(p,r\right)}=\sum_{k=1}^{n}H_{k}^{\left(p,r-1\right)},\text{
}\left(p\in\mathbb{Z},r\in\mathbb{N}\right),$ (8)

with $H_{n}^{\left(p,0\right)}=1/n^{p}$. These are a unified extension of
generalized harmonic numbers and hyperharmonic numbers:

$H_{n}^{\left(p,1\right)}=H_{n}^{\left(p\right)}\text{ and
}H_{n}^{\left(1,r\right)}=h_{n}^{(r)}.$

The main objective of this study is the evaluation of Euler sums of
generalized hyperharmonic numbers

$\zeta_{H^{\left(p,q\right)}}\left(r\right)=\sum\limits_{n=1}^{\infty}\dfrac{H_{n}^{\left(p,q\right)}}{n^{r}}.$

A step towards the solution of this problem is taken in [21]. However, the
recurrence used by the authors did not return a closed formula. Without an
available closed formula, they listed only the following few special cases

$\displaystyle\zeta_{H^{\left(p,2\right)}}\left(r\right)$
$\displaystyle=\zeta_{H^{\left(r-1\right)}}\left(p\right)-\zeta_{H^{\left(r\right)}}\left(p-1\right)+\zeta_{H^{\left(r\right)}}\left(p\right),$
$\displaystyle 2\zeta_{H^{\left(p,3\right)}}\left(r\right)$
$\displaystyle=2\zeta_{H^{\left(r\right)}}\left(p\right)+3\zeta_{H^{\left(r-1\right)}}\left(p\right)+\zeta_{H^{\left(r-2\right)}}\left(p\right)-3\zeta_{H^{\left(r\right)}}\left(p-1\right)$
$\displaystyle\quad+\zeta_{H^{\left(r\right)}}\left(p-2\right)-2\zeta_{H^{\left(r-1\right)}}\left(p-1\right).$

Later, Göral and Sertbaş [23] showed that the Euler sums of generalized
hyperharmonic numbers can be evaluated in terms of the Euler sums of
generalized harmonic numbers and special values of the Riemann zeta function.
However, their method does not determine the coefficients explicitly. This gap
is filled in this study. The following recurrence relation for
$H_{n}^{\left(p,q\right)}$ depending on the index $q$,

$\left(q-1\right)H_{n}^{\left(p,q\right)}=\left(n+q-1\right)H_{n}^{\left(p,q-1\right)}-H_{n}^{\left(p-1,q-1\right)}$

is obtained. Thanks to this recurrence relation, it is managed to obtain a
closed formula for $H_{n}^{\left(p,q\right)}$ in terms of
$H_{n}^{\left(p\right)}$ in Theorem 2. This enables the evaluation of Euler
sums of generalized hyperharmonic numbers in terms of the Euler sums of
generalized harmonic numbers as

$\zeta_{H^{\left(p,q+1\right)}}\left(r\right)=\frac{1}{q!}\sum_{m=0}^{q}\sum_{k=0}^{m}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q+1}{m+1}\binom{m}{k}\zeta_{H^{\left(p-k\right)}}\left(r+k-m\right).$

A demonstration of this formula is the following example:

$\displaystyle 2\zeta_{H^{\left(6,4+1\right)}}\left(6\right)$
$\displaystyle=-1925\zeta\left(11\right)+\left(175\pi^{2}-\frac{905}{4}-\frac{3937\pi^{8}}{544\,320}\right)\zeta\left(9\right)$
$\displaystyle+\left(\frac{245\pi^{2}}{12}+\frac{35\pi^{4}}{18}+\frac{31\pi^{10}}{46\,656}\right)\zeta\left(7\right)+\left(\frac{\pi^{4}}{4}+\frac{5\pi^{6}}{1134}+\frac{31\pi^{12}}{6123\,600}\right)\zeta\left(5\right)$
$\displaystyle-\frac{35}{12}\zeta^{2}\left(5\right)-\frac{1}{3}\zeta\left(3\right)\zeta\left(5\right)+\frac{\pi^{6}}{1134}\zeta\left(3\right)+\frac{\pi^{10}}{29\,160}+\frac{1406\pi^{12}}{638\,512\,875}.$

In addition, a counterpart of the reflection formula (4) is obtained in the
following form:

$\zeta_{H^{\left(p,q+1\right)}}\left(r\right)+\zeta_{H^{\left(r,q+1\right)}}\left(p\right).$

Section two completes with this formula which serves to calculate sums similar
to the foregoing example with less computational cost.

In the last section we further extend our results. In this direction we
establish new and more general identities for the series whose terms are
generalizations of harmonic numbers and reciprocal binomial coefficients. For
instance,

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,q\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$

is evaluated in terms of Riemann zeta values. This leads to several new
evaluation formulas for particular series involving generalized harmonic and
hyperharmonic numbers. We point out special cases of these formulas which
match with several known results in the literature.

##  2 Euler sums of generalized hyperharmonic numbers

In this section we present an evaluation formula for Euler sums
$\zeta_{H^{\left(p,q\right)}}\left(r\right)$ under certain conditions. To
state and prove our result we need some preliminaries.

Firstly, recall the polylogarithm defined by

$Li_{p}\left(t\right)=\sum_{k=1}^{\infty}\frac{t^{k}}{k^{p}},\text{
}(\left|t\right|\leq 1\text{ if }p>1,\text{ and }\left|t\right|<1\text{ if
}p\leq 1)\text{.}$

The generating function of the numbers $H_{n}^{\left(p,q\right)}$ in terms of
the polylogarithm is [21]

$\sum_{n=0}^{\infty}H_{n}^{\left(p,q\right)}t^{n}=\frac{Li_{p}\left(t\right)}{\left(1-t\right)^{q}},\text{
}\left|t\right|<1,\text{ }p,q\in\mathbb{Z}.$ (9)

Our first result presents the following reduction formula for
$H_{n}^{\left(p,q\right)}$.

###### Lemma 1

Let $p$ and $q$ be integers with $q\geq 1.$ Reduction relation for
$H_{n}^{\left(p,q\right)}$ in the index $q$ is

$\left(q-1\right)H_{n}^{\left(p,q\right)}=\left(n+q-1\right)H_{n}^{\left(p,q-1\right)}-H_{n}^{\left(p-1,q-1\right)}.$
(10)

Proof. We define the polynomial $H_{n}^{\left(p,q\right)}\left(z\right)$ as

$H_{n}^{\left(p,q\right)}\left(z\right)=\sum_{k=0}^{n}H_{k}^{\left(p,q\right)}z^{k}.$

Considering (9), we obtain the ordinary generating function of
$H_{n}^{\left(p,q\right)}\left(z\right)$ as

$\sum_{n=0}^{\infty}H_{n}^{\left(p,q\right)}\left(z\right)t^{n}=\frac{Li_{p}\left(zt\right)}{\left(1-t\right)\left(1-zt\right)^{q}}.$
(11)

From (11), it can be seen that

$z\frac{d}{dz}H_{n}^{\left(p,q\right)}\left(z\right)=H_{n}^{\left(p-1,q\right)}\left(z\right)+qzH_{n-1}^{\left(p,q+1\right)}\left(z\right).$
(12)

On the other hand, we utilize (8) twice to find that

$\displaystyle H_{n}^{\left(p,q\right)}$
$\displaystyle=\sum_{k=1}^{n}H_{k}^{\left(p,q-1\right)}=\sum_{k=1}^{n}\sum_{j=1}^{k}H_{j}^{\left(p,q-2\right)}$
$\displaystyle=\left(n+1\right)H_{n}^{\left(p,q-1\right)}-\sum_{j=1}^{n}jH_{j}^{\left(p,q-2\right)}$
$\displaystyle=\left(n+1\right)H_{n}^{\left(p,q-1\right)}-\left.\frac{d}{dz}H_{n}^{\left(p,q-2\right)}\left(z\right)\right|_{z=1},$

or equivalently

$\left.\frac{d}{dz}H_{n}^{\left(p,q-2\right)}\left(z\right)\right|_{z=1}=\left(n+1\right)H_{n}^{\left(p,q-1\right)}-H_{n}^{\left(p,q\right)}\overset{\text{(\ref{0})}}{=}nH_{n}^{\left(p,q-1\right)}-H_{n-1}^{\left(p,q\right)}.$
(13)

Therefore, (12) and (13) yield the desired formula.

The objective here is to express $H_{n}^{\left(p,q\right)}$ in terms of
$H_{n}^{\left(p\right)}.$ In [21] this relation is listed for at most $q=4$
due to the complexity of the process. However, the next result provides a
general solution to this problem where the numbers $H_{n}^{\left(p,q\right)}$
are expressed in terms of the numbers $H_{n}^{\left(p\right)}$ and
$\genfrac{[}{]}{0.0pt}{}{q}{j}_{r}$. Here $\genfrac{[}{]}{0.0pt}{}{q}{j}_{r}$
denotes the $r$-Stirling number of the first kind defined by the ”horizontal”
generating function [11, 12, 28]

$\left(x+r\right)\left(x+r+1\right)\cdots\left(x+r+q-1\right)=\sum\limits_{j=0}^{q}\genfrac{[}{]}{0.0pt}{}{q}{j}_{r}x^{j}.$
(14)

The essence of the theorem’s proof is based on the relationship between
$r$-Stirling numbers and symmetric polynomials. The $k$th elementary symmetric
polynomial $e_{k}\left(X_{1},\ldots,X_{q}\right)$ in variables
$X_{1},\ldots,X_{q}$ is defined by (see for example [27])

$\displaystyle e_{0}\left(X_{1},\ldots,X_{q}\right)$ $\displaystyle=1,$
$\displaystyle e_{k}\left(X_{1},\ldots,X_{q}\right)$
$\displaystyle=\sum_{1\leq j_{1}<j_{2}<\cdots<j_{k}\leq
q}\prod\limits_{i=1}^{k}X_{j_{i}},\text{ }1\leq k\leq q,$ $\displaystyle
e_{k}\left(X_{1},\ldots,X_{q}\right)$ $\displaystyle=0,\text{ }k>q,$

and possess the identity

$\prod\limits_{j=1}^{q}\left(x-X_{j}\right)=\sum\limits_{j=0}^{q}\left(-1\right)^{r-j}e_{q-j}\left(X_{1},\ldots,X_{q}\right)x^{j}.$
(15)

The comparison of (14) with (15) obviously leads to the following relationship
[28, Theorem 4.1]

$e_{q-j}\left(n+1,n+2,\ldots,n+q\right)=\genfrac{[}{]}{0.0pt}{}{q}{j}_{n+1}.$
(16)

###### Theorem 2

Let $p$ and $q$ be integers with $q\geq 0.$ Then,

$q!H_{n}^{\left(p,q+1\right)}=\sum_{k=0}^{q}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q}{k}_{n+1}H_{n}^{\left(p-k\right)}.$
(17)

Proof. We employ (10) on the right-hand side of

$\left(q-1\right)qH_{n}^{\left(p,q+1\right)}=\left(n+q\right)\left(q-1\right)H_{n}^{\left(p,q\right)}-\left(q-1\right)H_{n}^{\left(p-1,q\right)},$

and see that

$\displaystyle\left(q-1\right)qH_{n}^{\left(p,q+1\right)}$
$\displaystyle=H_{n}^{\left(p,q-1\right)}\left\\{\left(n+q\right)\left(n+q-1\right)\right\\}$
$\displaystyle\quad-
H_{n}^{\left(p-1,q-1\right)}\left\\{\left(n+q\right)+\left(n+q-1\right)\right\\}+H_{n}^{\left(p-2,q-1\right)}$
$\displaystyle=\sum_{k=0}^{2}\left(-1\right)^{k}e_{2-k}\left(n+q-1,n+q\right)H_{n}^{\left(p-k,q+1-2\right)}.$

These initial steps suggest that the following equality should hold:

$\displaystyle\left(q+1-r\right)\left(q+1-\left(r-1\right)\right)\cdots\left(q-1\right)qH_{n}^{\left(p,q+1\right)}$
(18)
$\displaystyle=\sum_{k=0}^{r}\left(-1\right)^{k}e_{r-k}\left(n+q-\left(r-1\right),n+q-\left(r-2\right),\ldots,n+q\right)H_{n}^{\left(p-k,q+1-r\right)}.$

To prove this by induction we show that it is also true for $r+1\leq q.$ We
multiply (18) by $\left(q-r\right)$ and then use (10). Hence we find that

$\displaystyle\left(q-r\right)\left(q+1-r\right)\left(q+1-\left(r-1\right)\right)\cdots\left(q-1\right)qH_{n}^{\left(p,q+1\right)}$
$\displaystyle\
=\sum_{k=0}^{r}\left(-1\right)^{k}e_{r-k}\left(n+q-\left(r-1\right),\ldots,n+q\right)\left(n+q-r\right)H_{n}^{\left(p-k,q-r\right)}$
$\displaystyle\quad+\sum_{k=0}^{r}\left(-1\right)^{k+1}e_{r-k}\left(n+q-\left(r-1\right),\ldots,n+q\right)H_{n}^{\left(p-k-1,q-r\right)}$
$\displaystyle\
=e_{r}\left(n+q-\left(r-1\right),\ldots,n+q\right)\left(n+q-r\right)H_{n}^{\left(p,q-r\right)}$
$\displaystyle\quad+\sum_{k=1}^{r}\left(-1\right)^{k}H_{n}^{\left(p-k,q-r\right)}\left\\{\left(n+q-r\right)e_{r-k}\left(n+q-\left(r-1\right),\ldots,n+q\right)\right.$
$\displaystyle\qquad\qquad\qquad\qquad\qquad\qquad\qquad\left.+e_{r+1-k}\left(n+q-\left(r-1\right),\ldots,n+q\right)\right\\}$
$\displaystyle\quad+\left(-1\right)^{r+1}e_{0}\left(n+q-r,\ldots,n+q\right)H_{n}^{\left(p-\left(r+1\right),q-r\right)}$
$\displaystyle\
=\sum_{k=0}^{r+1}\left(-1\right)^{k}e_{r+1-k}\left(n+q-r,\ldots,n+q\right)H_{n}^{\left(p-k,q+1-\left(r+1\right)\right)}.$

The case $r=q$ in (18) gives

$q!H_{n}^{\left(p,q+1\right)}=\sum_{k=0}^{q}\left(-1\right)^{k}e_{q-k}\left(n+1,\ldots,n+q\right)H_{n}^{\left(p-k\right)},$

which combines with (16) to give the statement.

Now, we are ready to state and prove our evaluation formula for
$\zeta_{H^{\left(p,q\right)}}\left(r\right).$ Thanks to this formula the
evaluation of Euler sums of generalized hyperharmonic numbers reduces to the
evaluation of Euler sums of generalized harmonic numbers.

###### Theorem 3

For $p,q\geq 1$ and $r>q+1,$ we have

$\zeta_{H^{\left(p,q+1\right)}}\left(r\right)=\frac{1}{q!}\sum_{m=0}^{q}\sum_{k=0}^{m}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q+1}{m+1}\binom{m}{k}\zeta_{H^{\left(p-k\right)}}\left(r+k-m\right).$

Proof. From (17) and the following identity [31, p. 1661]

$\genfrac{[}{]}{0.0pt}{}{n}{k}_{r+1}=\sum_{m=k}^{n}\genfrac{[}{]}{0.0pt}{}{n+1}{m+1}\binom{m}{k}r^{m-k},$

we have

$\displaystyle H_{n}^{\left(p,q+1\right)}$
$\displaystyle=\frac{1}{q!}\sum_{k=0}^{q}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q}{k}_{n+1}H_{n}^{\left(p-k\right)}$
$\displaystyle=\frac{1}{q!}\sum_{k=0}^{q}\sum_{m=k}^{q}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q+1}{m+1}\binom{m}{k}n^{m-k}H_{n}^{\left(p-k\right)}.$

Multiplying both sides with $n^{-r}$ and summing over $n$ from $1$ to
$\infty$, we deduce the desired result.

As mentioned introductory the sums
$\zeta_{H^{\left(p,q\right)}}\left(r\right)$ were listed up to $q=3$ in [21].
With the help of Theorem 3 these sums can be evaluated for further choices of
$q.$ For instance for $q=4$ one can obtain:

$\displaystyle\zeta_{H^{\left(p,4\right)}}\left(r\right)$
$\displaystyle=\zeta_{H^{\left(p\right)}}\left(r\right)+\frac{11}{6}\zeta_{H^{\left(p\right)}}\left(r-1\right)+\zeta_{H^{\left(p\right)}}\left(r-2\right)+\frac{1}{6}\zeta_{H^{\left(p\right)}}\left(r-3\right)$
$\displaystyle\quad-\frac{11}{6}\zeta_{H^{\left(p-1\right)}}\left(r\right)-2\zeta_{H^{\left(p-1\right)}}\left(r-1\right)-\frac{1}{2}\zeta_{H^{\left(p-1\right)}}\left(r-2\right)+\zeta_{H^{\left(p-2\right)}}\left(r\right)$
$\displaystyle\quad+\frac{1}{2}\zeta_{H^{\left(p-2\right)}}\left(r-1\right)-\frac{1}{6}\zeta_{H^{\left(p-3\right)}}\left(r\right).$

Hence, with the use of some values of
$\zeta_{H^{\left(p\right)}}\left(r\right)$ listed in forthcoming Remark 5, a
few concrete expressions of $\zeta_{H^{\left(p,4\right)}}\left(r\right)$ are:

$\displaystyle\bullet\text{ }\zeta_{H^{\left(1,4\right)}}\left(5\right)$
$\displaystyle=\frac{11}{2}\zeta\left(5\right)-\left(1-\frac{11}{36}\pi^{2}\right)\zeta\left(3\right)-\frac{1}{2}\left(\zeta\left(3\right)\right)^{2}-\frac{11}{216}\pi^{2}-\frac{\pi^{4}}{810}+\frac{\pi^{6}}{540},$
$\displaystyle\bullet\text{ }\zeta_{H^{\left(2,4\right)}}\left(5\right)$
$\displaystyle=-10\zeta\left(7\right)+\left(\frac{5}{6}\pi^{2}-\frac{21}{2}\right)\zeta\left(5\right)+\left(\frac{\pi^{4}}{45}+\frac{5}{6}\pi^{2}+\frac{5}{12}\right)\zeta\left(3\right)$
$\displaystyle\quad+\frac{11}{4}\left(\zeta\left(3\right)\right)^{2}+\frac{7\pi^{4}}{1080}-\frac{55\pi^{6}}{13608},$
$\displaystyle\bullet\text{ }\zeta_{H^{\left(3,4\right)}}\left(5\right)$
$\displaystyle=\zeta_{H^{\left(3\right)}}\left(5\right)+\frac{154}{3}\zeta\left(7\right)+\left(\frac{14}{3}-\frac{55}{12}\pi^{2}\right)\zeta\left(5\right)-2\left(\zeta\left(3\right)\right)^{2}$
$\displaystyle\quad-\left(\frac{7\pi^{2}}{18}+\frac{11\pi^{4}}{270}\right)\zeta\left(3\right)-\frac{\pi^{4}}{540}+\frac{\pi^{6}}{324},$
$\displaystyle\bullet\text{ }\zeta_{H^{\left(4,4\right)}}\left(5\right)$
$\displaystyle=-\frac{11}{6}\zeta_{H^{\left(3\right)}}\left(5\right)-\frac{125}{2}\zeta\left(9\right)+\left(\frac{35}{6}\pi^{2}-63\right)\zeta\left(7\right)$
$\displaystyle\quad+\left(\frac{35}{6}\pi^{2}+\frac{\pi^{4}}{18}\right)\zeta\left(5\right)+\frac{\pi^{4}}{30}\zeta\left(3\right)-\frac{\pi^{6}}{1944}+\frac{143\pi^{8}}{680400},$
$\displaystyle\bullet\text{ }\zeta_{H^{\left(5,4\right)}}\left(5\right)$
$\displaystyle=231\zeta\left(9\right)+\left(21-\frac{385}{18}\pi^{2}\right)\zeta\left(7\right)-\left(\frac{11}{60}\pi^{4}+\frac{23}{12}\pi^{2}\right)\zeta\left(5\right)$
$\displaystyle\quad+\frac{1}{2}\left(\zeta\left(5\right)\right)^{2}+\zeta\left(3\right)\zeta\left(5\right)-\frac{7\pi^{4}}{540}\zeta\left(3\right)-\frac{\pi^{8}}{8100}+\frac{\pi^{10}}{187110}.$

The following corollary gives the reflection formula for Euler sums of
generalized hyperharmonic numbers. Combined with (3), this corollary shows
that
$\zeta_{H^{\left(p,q+1\right)}}\left(r\right)+\zeta_{H^{\left(r,q+1\right)}}\left(p\right)$
can be written as a combination of Riemann zeta values. In this way,
particular Euler sums of type $\zeta_{H^{\left(p,q\right)}}\left(p\right)$ can
be evaluated with less computation.

###### Corollary 4

Let $p>q+1,$ $r>q+1$ and $p+r$ be even. Then

$\displaystyle\zeta_{H^{\left(p,q+1\right)}}\left(r\right)+\zeta_{H^{\left(r,q+1\right)}}\left(p\right)$
$\displaystyle=\zeta\left(p+r\right)+\frac{2}{q!}\sum_{\begin{subarray}{c}m=0\\\
m\text{
odd}\end{subarray}}^{q}\sum_{k=0}^{m}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q+1}{m+1}\binom{m}{k}\zeta_{H^{\left(p-k\right)}}\left(r+k-m\right)$
$\displaystyle+\frac{1}{q!}\sum_{m=0}^{q}\sum_{k=0}^{m}\left(-1\right)^{m+k}\genfrac{[}{]}{0.0pt}{}{q+1}{m+1}\binom{m}{k}\zeta\left(p-k\right)\zeta\left(r+k-m\right).$

Proof. Let $\left(p+r\right)$ be even. It is obvious from Theorem 3 that

$\displaystyle\zeta_{H^{\left(p,q+1\right)}}\left(r\right)+\zeta_{H^{\left(r,q+1\right)}}\left(p\right)$
$\displaystyle=\frac{1}{q!}\sum_{m=0}^{q}\sum_{k=0}^{m}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q+1}{m+1}\binom{m}{k}\left\\{\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p-k\right)}}{n^{r+k-m}}+\left(-1\right)^{m}\sum_{n=1}^{\infty}\frac{H_{n}^{\left(r+k-m\right)}}{n^{p-k}}\right\\}.$

We write the right-hand side as

$\displaystyle\sum_{\begin{subarray}{c}m=0\\\ m\text{
odd}\end{subarray}}^{q}\sum_{k=0}^{m}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q+1}{m+1}\binom{m}{k}\left\\{\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p-k\right)}}{n^{r+k-m}}-\sum_{n=1}^{\infty}\frac{H_{n}^{\left(r+k-m\right)}}{n^{p-k}}\right\\}$
$\displaystyle+\sum_{0\leq m\leq
q/2}\sum_{k=0}^{2m}\left(-1\right)^{k}\genfrac{[}{]}{0.0pt}{}{q+1}{2m+1}\binom{2m}{k}\left\\{\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p-k\right)}}{n^{r+k-2m}}+\sum_{n=1}^{\infty}\frac{H_{n}^{\left(r+k-2m\right)}}{n^{p-k}}\right\\}.$

By the reflection formula (4) we have

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p-k\right)}}{n^{r+k-2m}}+\sum_{n=1}^{\infty}\frac{H_{n}^{\left(r+k-2m\right)}}{n^{p-k}}=\zeta\left(p+r-2m\right)+\zeta\left(p-k\right)\zeta\left(r+k-2m\right).$

Moreover, for odd $m,$ it can be seen from (3) that

$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p-k\right)}}{n^{r+k-m}}-\sum_{n=1}^{\infty}\frac{H_{n}^{\left(r+k-m\right)}}{n^{p-k}}$
$\displaystyle=2\zeta_{H^{\left(p-k\right)}}\left(r+k-m\right)-\zeta\left(p+r-m\right)$
$\displaystyle\quad-\zeta\left(p-k\right)\zeta\left(r+k-m\right).$

Hence, we obtain the desired equation.

###### Remark 5

For interested readers we would like to list some values of
$\zeta_{H^{\left(p\right)}}\left(r\right)$, used in the evaluations of
$\zeta_{H^{\left(p,4\right)}}\left(5\right),$ $1\leq p\leq 5,$ and
$\zeta_{H^{\left(6,5\right)}}\left(6\right).$ These are calculated with the
help of (1), (3) and (4).

$\bullet\text{ }\zeta_{H^{\left(1\right)}}\left(2\right)=2\zeta\left(3\right),$ | $\bullet\text{ }\zeta_{H^{\left(3\right)}}\left(6\right)=\frac{85}{2}\zeta\left(9\right)-\frac{7\pi^{2}}{2}\zeta\left(7\right)-\frac{\pi^{4}}{15}\zeta\left(5\right),$
---|---
$\bullet\text{ }\zeta_{H^{\left(1\right)}}\left(3\right)=\frac{\pi^{4}}{72},$ | $\bullet\text{ }\zeta_{H^{\left(4\right)}}\left(2\right)=-\zeta^{2}\left(3\right)+\frac{37\pi^{6}}{11340},$
$\bullet\text{ }\zeta_{H^{\left(1\right)}}\left(4\right)=3\zeta\left(5\right)-\frac{\pi^{2}}{6}\zeta\left(3\right),$ | $\bullet\text{ }\zeta_{H^{\left(4\right)}}\left(3\right)=-17\zeta\left(7\right)+\frac{5\pi^{2}}{3}\zeta\left(5\right)+\frac{\pi^{4}}{90}\zeta\left(3\right),$
$\bullet\text{ }\zeta_{H^{\left(1\right)}}\left(5\right)=-\frac{1}{2}\zeta^{2}\left(3\right)+\frac{\pi^{6}}{540},$ | $\bullet\text{ }\zeta_{H^{\left(4\right)}}\left(4\right)=\frac{13\pi^{8}}{113\,400},$
$\bullet\text{ }\zeta_{H^{\left(2\right)}}\left(2\right)=\frac{7\pi^{4}}{360},$ | $\bullet\text{ }\zeta_{H^{\left(4\right)}}\left(5\right)=-\frac{125}{2}\zeta\left(9\right)+\frac{35\pi^{2}}{6}\zeta\left(7\right)+\frac{\pi^{4}}{18}\zeta\left(5\right),$
$\bullet\text{ }\zeta_{H^{\left(2\right)}}\left(3\right)=-\frac{9}{2}\zeta\left(5\right)+\frac{\pi^{2}}{2}\zeta\left(3\right),$ | $\bullet\text{ }\zeta_{H^{\left(5\right)}}\left(2\right)=11\zeta\left(7\right)-\frac{2\pi^{2}}{3}\zeta\left(5\right)-\frac{\pi^{4}}{45}\zeta\left(3\right),$
$\bullet\text{ }\zeta_{H^{\left(2\right)}}\left(4\right)=\zeta^{2}\left(3\right)-\frac{\pi^{6}}{2835},$ | $\bullet\text{ }\zeta_{H^{\left(5\right)}}\left(4\right)=\frac{127}{2}\zeta\left(9\right)-\frac{35\pi^{2}}{6}\zeta\left(7\right)-\frac{2\pi^{4}}{45}\zeta\left(5\right),$
$\bullet\text{ }\zeta_{H^{\left(2\right)}}\left(5\right)=-10\zeta\left(7\right)+\frac{5\pi^{2}}{6}\zeta\left(5\right)$ | $\bullet\text{ }\zeta_{H^{\left(5\right)}}\left(5\right)=\frac{1}{2}\zeta^{2}\left(5\right)+\frac{\pi^{10}}{187110},$
$\qquad\qquad\qquad+\frac{\pi^{4}}{45}\zeta\left(3\right),$ | $\bullet\text{ }\zeta_{H^{\left(5\right)}}\left(6\right)=\frac{463}{2}\zeta\left(11\right)-21\pi^{2}\zeta\left(9\right)-\frac{7}{30}\pi^{4}\zeta\left(7\right),$
$\bullet\text{ }\zeta_{H^{\left(3\right)}}\left(2\right)=\frac{11}{2}\zeta\left(5\right)-\frac{\pi^{2}}{3}\zeta\left(3\right),$ | $\bullet\text{ }\zeta_{H^{\left(6\right)}}\left(3\right)=\frac{7\pi^{2}}{2}\zeta\left(7\right)-\frac{83}{2}\zeta\left(9\right)+\frac{\pi^{4}}{15}\zeta\left(5\right)$
$\bullet\text{ }\zeta_{H^{\left(3\right)}}\left(3\right)=\frac{1}{2}\zeta^{2}\left(3\right)+\frac{\pi^{6}}{1890},$ | $\qquad\qquad\qquad+\frac{\pi^{6}}{945}\zeta\left(3\right),$
$\bullet\text{ }\zeta_{H^{\left(3\right)}}\left(4\right)=18\zeta\left(7\right)-\frac{5\pi^{2}}{3}\zeta\left(5\right),$ | $\bullet\text{ }\zeta_{H^{\left(6\right)}}\left(5\right)=21\pi^{2}\zeta\left(9\right)-\frac{461}{2}\zeta\left(11\right)+\frac{7\pi^{4}}{30}\zeta\left(7\right)$
| $\qquad\qquad\qquad+\frac{\pi^{6}}{945}\zeta\left(5\right).$

##  3 Series involving harmonic numbers and reciprocal binomial coefficients

In this section, we introduce evaluation formulas for some series involving
the harmonic numbers and their generalizations.

###### Theorem 6

Let $p\geq 1$ and $q,$ $l\geq 0$ be integers with $l\geq q$. For $m\geq 1$,

$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,q\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
(19) $\displaystyle\
=\sum_{j=0}^{l-q}\binom{l-q}{j}\left\\{\frac{\left(-1\right)^{j+p-1}H_{m+j}}{\left(m+j\right)^{p}}+\sum_{k=1}^{p-1}\frac{\left(-1\right)^{j+k-1}}{\left(m+j\right)^{k}}\zeta\left(p+1-k\right)\right\\}$

and

$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,q\right)}}{n\binom{n+l}{l}}$
(20) $\displaystyle\
=\zeta\left(p+1\right)-\sum_{j=1}^{l-q}\binom{l-q}{j}\left\\{\frac{\left(-1\right)^{j+p}H_{j}}{j^{p}}+\sum_{k=1}^{p-1}\frac{\left(-1\right)^{j+k}}{j^{k}}\zeta\left(p+1-k\right)\right\\}.$

Proof. Using the formula (see [25, p.909])

${\displaystyle\int\limits_{0}^{1}}t^{n+m-1}\left(1-t\right)^{l}dt=\frac{1}{\left(n+m\right)\binom{n+m+l}{l}},$

we can write

$\frac{H_{n}^{\left(p,q\right)}}{\left(n+m\right)\binom{n+m+l}{l}}={\displaystyle\int\limits_{0}^{1}}H_{n}^{\left(p,q\right)}t^{n+m-1}\left(1-t\right)^{l}dt.$

With the help of (9), we get

$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,q\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
$\displaystyle={\displaystyle\int\limits_{0}^{1}}t^{m-1}\left(1-t\right)^{l-q}Li_{p}\left(t\right)dt$
$\displaystyle=\sum_{j=0}^{l-q}\binom{l-q}{j}\left(-1\right)^{j}\sum_{n=1}^{\infty}\frac{1}{n^{p}\left(n+m+j\right)}.$

Then (19) follows from (7). If $m=0,$ then we have

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,q\right)}}{n\binom{n+l}{l}}=\zeta\left(p+1\right)+\sum_{j=1}^{l-q}\binom{l-q}{j}\left(-1\right)^{j}\sum_{n=1}^{\infty}\frac{1}{n^{p}\left(n+j\right)},$

which is equivalent to (20).

Now, we deal with some special cases of Theorem 6. Setting $q=l$ gives

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,q\right)}}{\left(n+m\right)\binom{n+m+q}{q}}=\frac{\left(-1\right)^{p-1}H_{m}}{m^{p}}+\sum_{i=1}^{p-1}\frac{\left(-1\right)^{i-1}}{m^{i}}\zeta\left(p+1-i\right)$
(21)

and

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,q\right)}}{n\binom{n+q}{q}}=\zeta\left(p+1\right).$
(22)

Note that the variable $q$ does not appear in the right-hand sides and all
these series converge very slowly.

For $p=1,$ (21) and (22) give Proposition 5 and Proposition 6 in [20]

$\sum_{n=1}^{\infty}\frac{h_{n}^{\left(q\right)}}{n\binom{n+q}{q}}=\frac{1}{6}\pi^{2}\text{
\ \ and \ \
}\sum_{n=1}^{\infty}\frac{h_{n-1}^{\left(q\right)}}{n\binom{n+q}{q}}=1,$

respectively. (21) also yields [44, Eq. (2.30)] for $q=1$. Additionally, when
$p=1$ in Theorem 6, we reach that

$\sum_{n=1}^{\infty}\frac{h_{n}^{\left(q\right)}}{\left(n+m\right)\binom{n+m+l}{l}}=\sum_{j=0}^{l-q}\binom{l-q}{j}\left(-1\right)^{j}\frac{H_{m+j}}{m+j}$

and

$\sum_{n=1}^{\infty}\frac{h_{n}^{\left(q\right)}}{n\binom{n+l}{l}}=\frac{1}{6}\pi^{2}+\sum_{j=1}^{l-q}\left(-1\right)^{j}\binom{l-q}{j}\frac{H_{j}}{j}.$

Now employing [15, Eq.(18)]

$\sum_{k=0}^{m}\left(-1\right)^{k}\binom{m}{k}\frac{H_{n+k}}{n+k}=\frac{H_{n+m}-H_{m}}{n\binom{n+m}{m}}$

and [9, Eq. (9.4b)]

$\sum_{j=1}^{n}\left(-1\right)^{j+1}\binom{n}{j}\frac{H_{j}}{j}=H_{n}^{\left(2\right)}$

gives the following closed forms for series involving hyperharmonic numbers
with reciprocal binomial coefficients.

###### Corollary 7

Let $q,l\geq 0$ be integers with $l\geq q$. For all integers $m\geq 1$

$\sum_{n=1}^{\infty}\frac{h_{n}^{\left(q\right)}}{\left(n+m\right)\binom{n+m+l}{l}}=\frac{H_{m+l-q}-H_{l-q}}{m\binom{n+m+l-q}{l-q}}$

and

$\sum_{n=1}^{\infty}\frac{h_{n}^{\left(q\right)}}{n\binom{n+l}{l}}=\frac{1}{6}\pi^{2}-H_{l-q}^{\left(2\right)}.$

For $p=q=1,$ Theorem 6 gives [39, Eq. (2.31)]

$\sum_{n=1}^{\infty}\frac{H_{n}}{\left(n+m\right)\binom{n+m+q}{q}}=\frac{H_{n+q-1}-H_{q-1}}{n\binom{n+q-1}{q-1}}$

and [38]

$\sum_{n=1}^{\infty}\frac{H_{n}}{n\binom{n+l}{l}}=\frac{1}{6}\pi^{2}-H_{l-1}^{\left(2\right)}.$

For $q=1,$ (20) becomes

$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{n\binom{n+l}{l}}$
(23) $\displaystyle\
=\zeta\left(p+1\right)-\sum_{j=1}^{l-1}\binom{l-1}{j}\left\\{\frac{\left(-1\right)^{j+p}H_{j}}{j^{p}}+\sum_{k=1}^{p-1}\frac{\left(-1\right)^{j+k}}{j^{k}}\zeta\left(p+1-k\right)\right\\},$

which is also given by Sofo [33, Theorem 2.2] in a slightly different form.

Setting $q=1$ in (19) yields the following corollary involving generalized
harmonic numbers.

###### Corollary 8

For all integers $m,p,l\geq 1,$

$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
(24)
$\displaystyle=\sum_{j=0}^{l-1}\binom{l-1}{j}\left(-1\right)^{j}\left\\{\frac{\left(-1\right)^{p-1}H_{m+j}}{\left(m+j\right)^{p}}+\sum_{k=1}^{p-1}\frac{\left(-1\right)^{k-1}}{\left(m+j\right)^{k}}\zeta\left(p+1-k\right)\right\\}.$

The following particular cases can be deduced setting $p=2$ and $p=3$:

$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(2\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
$\displaystyle=\frac{\pi^{2}}{6m\binom{m+l-1}{l-1}}-\sum_{j=0}^{l-1}\left(-1\right)^{j}\binom{l-1}{j}\frac{H_{m+j}}{\left(m+j\right)^{2}},$
(25)
$\displaystyle\sum_{n=1}^{\infty}\frac{H_{n}^{\left(3\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
$\displaystyle=\sum_{j=0}^{l-1}\left(-1\right)^{j}\binom{l-1}{j}\frac{H_{m+j}}{\left(m+j\right)^{3}}$
$\displaystyle+\frac{1}{m\binom{m+l-1}{l-1}}\left\\{\zeta\left(3\right)-\frac{\pi^{2}}{6}\left(H_{m+l-1}-H_{m-1}\right)\right\\}.$

It is worth noting that the special case with choices $m=6$ and $l=3$ in (25)
is recorded in [35, Remark 1] despite a misprint. The case is below

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(2\right)}}{\left(n+6\right)\binom{n+9}{3}}=\frac{1}{168}\zeta\left(2\right)-\frac{37073}{7902720}.$

For our final results, we deal with the special case $q=2$ of Theorem 6. By
aid of (17),
$\genfrac{[}{]}{0.0pt}{}{q}{0}_{r}=r\left(r+1\right)\cdots\left(r+q-1\right)$
and $\genfrac{[}{]}{0.0pt}{}{q}{q}_{r}=1,$ we have

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{\binom{n+l}{l}}=\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,2\right)}}{n\binom{n+l}{l}}+\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p-1\right)}}{n\binom{n+l}{l}}-\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{n\binom{n+l}{l}},$

where $l$ is any integer greater than $1$. From (20), (23) and some
arrangements we obtain

$\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{\binom{n+l}{l}}=\zeta\left(p\right)+\sum_{j=1}^{l-1}\left(-1\right)^{j}\left\\{\binom{l-1}{j}\mu\left(p-1,j\right)-\binom{l-2}{j-1}\mu\left(p,j\right)\right\\},$

where $\mu\left(p,j\right)$ is given in (7). A slightly different form of the
equation above is given in [33, Theorem 2.1]. Similarly

$\displaystyle\sum_{n=1}^{\infty}\frac{nH_{n}^{\left(p\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
$\displaystyle=\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p,2\right)}}{\left(n+m\right)\binom{n+m+l}{l}}+\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p-1\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
$\displaystyle-\sum_{n=1}^{\infty}\frac{H_{n}^{\left(p\right)}}{\left(n+m\right)\binom{n+m+l}{l}}.$

Then exploiting (20) and (24) in the last equation yields the following
corollary.

###### Corollary 9

Let $l>1$ be an integer. Then

$\displaystyle\sum_{n=1}^{\infty}\frac{nH_{n}^{\left(p\right)}}{\left(n+m\right)\binom{n+m+l}{l}}$
$\displaystyle\
=\mu\left(p-1,m\right)+\sum_{j=1}^{l-1}\left(-1\right)^{j}\left\\{\binom{l-1}{j}\mu\left(p-1,m+j\right)-\binom{l-2}{j-1}\mu\left(p,m+j\right)\right\\},$

where $\mu\left(p,j\right)$ is given in (7).

For $l=2$ this formula can be read as

$\displaystyle\sum_{n=1}^{\infty}\frac{nH_{n}^{\left(p\right)}}{\left(n+m\right)\left(n+m+1\right)\left(n+m+2\right)}$
$\displaystyle=\left(-1\right)^{p+1}\frac{m+2}{2\left(m+1\right)^{p}}H_{m+1}+\left(-1\right)^{p}\frac{H_{m}}{2m^{p-1}}$
$\displaystyle\quad+\sum_{k=1}^{p-2}\left(-1\right)^{k-1}\left\\{\frac{1}{2m^{k}}-\frac{1}{2\left(m+1\right)^{k}}\right\\}\zeta\left(p-k\right)+\sum_{k=1}^{p-1}\frac{\left(-1\right)^{k-1}}{2\left(m+1\right)^{k}}\zeta\left(p+1-k\right).$

The first few cases of this formula are listed below:

$\displaystyle\sum_{n=1}^{\infty}\frac{nH_{n}}{\left(n+m\right)\left(n+m+1\right)\left(n+m+2\right)}$
$\displaystyle=\frac{1}{2\left(m+1\right)}\left(H_{m+1}+1\right),$
$\displaystyle\sum_{n=1}^{\infty}\frac{nH_{n}^{\left(2\right)}}{\left(n+m\right)\left(n+m+1\right)\left(n+m+2\right)}$
$\displaystyle=\frac{\pi^{2}}{12\left(m+1\right)}-\frac{\left(m+2\right)}{2\left(m+1\right)^{2}}H_{m+1}+\frac{1}{2m}H_{m},$
$\displaystyle\sum_{n=1}^{\infty}\frac{nH_{n}^{\left(3\right)}}{\left(n+m\right)\left(n+m+1\right)\left(n+m+2\right)}$
$\displaystyle=\frac{\zeta\left(3\right)}{2\left(m+1\right)}+\frac{\pi^{2}}{12m\left(m+1\right)^{2}}$
$\displaystyle\ \ \
+\frac{m+2}{2\left(m+1\right)^{3}}H_{m+1}-\frac{1}{2m^{2}}H_{m}.$

## References

  * [1] Adamchik V. On Stirling numbers and Euler sums, Journal of Computational and Applied Mathematics 1997; 79 (1): 119–130.
  * [2] Alzer H. and Choi J. Four parametric linear Euler sums. J. Math. Anal. Appl., 484 (2020), 123661.
  * [3] Bailey DH, Borwein D, Borwein JM. On Eulerian log-gamma integrals and Tornheim–Witten zeta functions. The Ramanujan Journal 2015; 36 (1-2): 43–68.
  * [4] Basu A. On the evaluation of Tornheim sums and allied double sums. The Ramanujan Journal 2011; 26 (2): 193–207.
  * [5] Berndt BC. Ramanujan’s Notebooks. Part I, New York, USA: Springer-Verlag, 1985.
  * [6] Borwein D, Borwein JM, Girgensohn R. Explicit evaluation of Euler sums. Proceedings of the Edinburgh Mathematical Society 1995; 38 (2): 277–294.
  * [7] Borwein JM, Girgensohn R. Evaluation of triple Euler sums. The Electronic Journal of Combinatorics 1996; 3: #R23.
  * [8] Boyadzhiev KN. Consecutive evaluation of Euler sums. International Journal of Mathematics and Mathematical Sciences 2002; 29 (9): 555–561.
  * [9] Boyadzhiev KN. Notes on the Binomial Transform, Theory and Table, World Scientific, 2018.
  * [10] Broadhurst DJ. On the enumeration of irreducible $k$-fold Euler sums and their roles in knot theory and field theory, 1996 arXiv:hepth/9604128.
  * [11] Broder AZ. The $r$-Stirling numbers, Discrete Mathematics 1984; 49: 241–259.
  * [12] L. Carlitz, Weighted Stirling numbers of the first and second kind-I, II, Fibonacci Quart. 18 (1980), 147-162, 242-257.
  * [13] Choi J. Summation formulas involving binomial coefficients, harmonic numbers, and generalized harmonic numbers. Abstract and Applied Analysis. Vol. 2014. Hindawi, 2014.
  * [14] J. Choi, H. Srivastava, Explicit Evaluation of Euler and Related Sums. Ramanujan J., 10 (2005), 51–70.
  * [15] Chu W. Summation formulae involving harmonic numbers. Filomat 2012; 26 (1): 143–152.
  * [16] Coffey MW. On some log-cosine integrals related to $\zeta\left(3\right)$, $\zeta\left(4\right)$ and $\zeta\left(6\right)$. Journal of Computational and Applied Mathematics 2003; 159 (2): 205–215.
  * [17] Coffey MW. On one-dimensional digamma and polygamma series related to the evaluation of Feynman diagrams. Journal of Computational and Applied Mathematics 2005; 183: 84–100.
  * [18] Conway JH, Guy RK. The Book of Numbers. New York, USA: Springer-Verlag, 1996.
  * [19] Crandall RE, BuhlerJP. On the evaluation of Euler sums. Experimental Mathematics 1994; 3: 275–285.
  * [20] Dil A, Boyadzhiev KN. Euler sums of hyperharmonic numbers. Journal of Number Theory 2015; 147: 490–498.
  * [21] Dil A, Mező I, Cenkci M. Evaluation of Euler-like sums via Hurwitz zeta values. Turkish Journal of Mathematics 2017; 41: 1640–1655.
  * [22] Euler L. Meditationes circa singulare serierum genus. Novi Commentarii academiae scientiarum Petropolitanae 1776; 20: 140–186.
  * [23] Göral H, Sertbaş DC. Euler sums and non-integerness of harmonic type sums. Hacettepe Journal of Mathematics and Statistics 2020; 49 (2): 586–598.
  * [24] Huard JG, Williams KS, Nan-Yue Z. On Tornheim’s double series, Acta Arithmetica 1996; 75 (2): 105–117.
  * [25] Jeffrey A, Zwillinger D. (Eds.) Gradshteyn and Ryzhik’s Table of Integrals, Series, and Products. $7$th edition, Burlington, USA: Academic Press, 2007.
  * [26] Kamano K. Dirichlet series associated with hyperharmonic numbers. Memoirs of the Osaka Institute of Technology, Series A 2011; 56: 11–15.
  * [27] Macdonald IG. Symmetric Functions and Hall Polynomials. Oxford Mathematical Monographs, New York, USA: Oxford University Press, 1995.
  * [28] Merris R. The p-Stirling numbers. Turkish Journal of Mathematics 2000; 24: 379–399.
  * [29] Mező I, Dil A. Hyperharmonic series involving Hurwitz zeta function. Journal of Number Theory 2010; 130 (2): 360–369.
  * [30] Nielsen N. Handbuch der Theorie der Gammafunktion. New York, USA: Reprinted by Chelsea Publishing Company, 1965.
  * [31] Nyul G, Rácz G. The r-Lah numbers, DiscreteMath. 338 (2015), 1660-1666.
  * [32] Sofo A. Harmonic sums and integral representations. Journal of Applied Analysis 2010; 16 (2): 265–277.
  * [33] Sofo A. Harmonic number sums in higher powers. Journal of Mathematical Analysis 2011; 2 (2): 15–22.
  * [34] Sofo A. Summation formula involving harmonic numbers. Analysis Mathematica 2011; 37 (1): 51–64.
  * [35] Sofo, A. Shifted harmonic sums of order two. Commun. Korean Math. Soc. 2014; 29 (2): 239–255.
  * [36] Sofo A. Polylogarithmic connections with Euler sums. Sarajevo Journal of Mathematics 2016; 12 (24): 17–32.
  * [37] Sofo A. General order Euler sums with multiple argument. Journal of Number Theory 2018; 189: 255–271.
  * [38] Sofo A., Cvijović D. Extensions of Euler harmonic sums. Applicable Analysis and Discrete Mathematics 2012: 6: 317-328.
  * [39] Sofo, A., Srivastava, H.M. A family of shifted harmonic sums. Ramanujan J. 2015; 37: 89–108.
  * [40] Srivastava HM, Choi J, Series associated with the zeta and related functions, Kluwer Academic Publishers, Dordrecht-Boston-London, 2001.
  * [41] Xu C. Euler sums of generalized hyperharmonic numbers. Journal of the Korean Mathematical Society 2018; 55 (5): 1207–1220.
  * [42] Xu C. Computation and theory of Euler sums of generalized hyperharmonic numbers. Comptes Rendus Mathematique 2018; 356 (3): 243–252.
  * [43] Xu C, Yan Y, Shi Z. Euler sums and integrals of polylogarithm functions, Journal of Number Theory 2016; 165: 84–108.
  * [44] Xu C, Zhang M, Zhu W. Some evaluation of harmonic number sums. Integral Transforms and Special Functions 2016; 27 (12): 937–955.
  * [45] Yang J, Wang Y. Summation formulae in relation to Euler sums. Integral Transforms and Special Functions 2017; 28 (5): 336–349.
  * [46] P. T. Young, Global series for zeta functions. Fibonacci Quart., 57 (2019), 154-169.
  * [47] P. T. Young, Absolute series for higher Euler constants. J. Number Theory, 220 (2021), 19-33.
  * [48] Wang W, Yanhong L. Euler sums and Stirling sums. Journal of Number Theory 2018; 185: 160–193.
  * [49] Zhou X, Cai T, Bradley D. Signed $q$-analogs of Tornheim’s double series. Proceedings of the American Mathematical Society 2008; 136 (8): 2689–2698.