# Aerodynamic effects and performance improvements of running
in drafting formations

Lukas Schickhofer  KTM E-TECHNOLOGIES, Salzburg, AT-5081, Austria
<EMAIL_ADDRESS>Henry Hanson  ADIDAS, Herzogenaurach, DE-91074,
Germany

###### Abstract

Drafting as a process to reduce drag and to benefit from the presence of other
competitors is applied in various sports with several recent examples of
competitive running in formations. In this study, the aerodynamics of a
realistic model of a female runner is calculated by computational fluid
dynamics (CFD) simulations at four running speeds of
$15\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$,
$18\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$,
$21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$, and
$36\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$. Aerodynamic power
fractions of the total energy expenditure are found to be in the range of
$2.6\%$-$8.5\%$. Additionally, four exemplary formations are analysed with
respect to their drafting potential and resulting drag values are compared for
the main runner and her pacers. The best of the formations achieves a total
drag reduction on the main runner of $75.6\%$. Moreover, there are large
variations in the drag reduction between the considered formations of up to
$42\%$ with respect to the baseline single-runner case. We conclude that major
drag reduction of more than $70\%$ can already be achieved with fairly simple
formations, while certain factors, such as runners on the sides, can have a
detrimental effect on drag reduction due to local acceleration of the passing
flow. Using an empirical model for mechanical power output during running,
gains of metabolic power and performance predictions are evaluated for all
considered formations. Improvements in running economy are up to $3.5\%$ for
the best formation, leading to velocity gains of $2.3\%$. This translates to
$154\text{\,}\mathrm{s}$ ($\approx$2.6\text{\,}\mathrm{min}$$) saved over a
marathon distance. Consequently, direct conclusions are drawn from the
obtained data for ideal drafting of long-distance running in highly packed
formations.

###### keywords:

Running aerodynamics, Drafting, Drag reduction, Metabolic power, Computational
fluid dynamics

###### MSC:

[2010] 00-01, 99-00

††journal: Journal of Biomechanics

##  1 Introduction

In most sports involving the race against other participants for time, drag
plays a crucial role. It is defined as the aerodynamic force acting in an
axial direction on the athlete as he or she moves against the surrounding
fluid, such as water or air. As a result, it has become common practice both
in sports and motorsports to benefit from the presence of other athletes or
vehicles to alleviate the resistance of the fluid by drafting in the
associated wake region (i.e. _slipstreaming_). In the wake there is not just
drastically lower fluid velocity, but also a negative pressure coefficient
leading to suction, which additionally benefits the trailing competitor.

While in motorsports with its high speeds and considerable research efforts,
drafting has long been investigated (e.g. Katz [1995], Romberg et al. [1971]),
it is less researched in other sports. However, Rundell [1996] investigated
drafting during speed skating and showed large positive effects on metabolic
activity, heart rate, and lactate response. Also in swimming, drafting was
associated with considerable energy saving and the optimal drafting distance
was established at $0$-$0.5\text{\,}\mathrm{m}$ behind the leading swimmer,
causing power savings of up to $31\%$ [Bassett Jr et al., 1991, Chatard and
Wilson, 2003]. Furthermore, cycling is a sport of large benefits through
aerodynamic drafting and riding in formations, such as tightly packed pelotons
[Edwards and Byrnes, 2007, Blocken et al., 2018, Malizia and Blocken, 2020].
Apart from testing in wind tunnels, typically applying similarity principles
through scaled models at similar Reynolds numbers, computational fluid
dynamics (CFD) has become a widely used method to accurately assess drag
distributions across athletes. Blocken et al. [2013] examined the effect of
two drafting cyclists at various upper-body positions through CFD and found
drag reductions of up to $27.1\%$ for the trailing cyclist, while also the
leading cyclist's drag was decreased by up to $2.6\%$. Moreover, drag
reductions of $42\%$ to $48\%$ have been measured for drafting during cycling
in the velodrome [Fitton et al., 2018]. Fuss [2018] recently suggested an
analytical model for slipstreaming in gravity-powered sports such as ski cross
and found significant advantages in gained velocity and glide distance with
respect to a leading skier, especially in a tucked body position.

Although running happens typically at fairly low speeds, drafting shows also
here a measurable impact on aerodynamic forces. In a major study on the
subject, Pugh [1970] demonstrated drag savings of up to $80\%$ when running at
a middle-distance speed of $4.46\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$
against a slight head wind. More recently, particularly in light of attempts
to break the $2$-hour barrier of the marathon distance (e.g. _INEOS 1:59
Challenge_ , _Nike Breaking2_), increasing research efforts were directed to
minimize drag for a main runner in formations. Hoogkamer et al. [2018] used a
reduced-order model to establish a sustainable velocity of
$5.93\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$ using cooperative drafting
in a four-runner team. Polidori et al. [2020] applied CFD to compute the drag
and power savings of Kenenisa Bekele, while running the Berlin marathon and
using cooperative drafting. It was shown that up to $57.3\%$ of aerodynamic
power were gained in the ideal case of running directly behind a front pacer
at a distance of $1.3\text{\,}\mathrm{m}$ when compared to the case of running
alone.

The current study aims at shedding more light on the aerodynamic effects of
running alone and in formations by using a realistic model of a female runner
and a validated CFD methodology. Four different formations are chosen to
reflect various running scenarios and their implications on the drag on the
main runner as well as her pacers are discussed. The effect of the resulting
drag reduction on the metabolic rate and running economy is computed by using
a mechanical power model of running. As a result, possible performance
predictions are made in terms of improvements of running speed and time.
Finally, by applying the conclusions from this study, direct implications for
an optimized drafting strategy are achieved.

##  2 Methods

###  2.1 Numerical method

For the computation of the flow field variables and the resulting aerodynamic
forces on the runner, computational fluid dynamics (CFD) using the finite-
volume method is applied. The Navier-Stokes equations describing the
conservation of momentum are discretized in its Reynolds-averaged (RANS) form.

Turbulence models have an important effect on aerodynamic force computations,
mainly because of their modelling of momentum transfer near surfaces and in
boundary layers. This in turn influences the locations of flow separation,
which is crucial for correct drag predictions. Here, turbulence is treated
with the $k$-$\omega$ model in its shear-stress-transport (SST) formulation
according to Menter [1994]. The model has been repeatedly shown to be suitable
for external aerodynamics and to give superior performance to the Spalart-
Allmaras and the $k$-$\epsilon$ models in flow scenarios with adverse pressure
gradients and free shear layers, which do occur in our setup [Menter, 1992,
Hellsten, 1998]. Defraeye et al. [2010] tested the performance of various
turbulence treatments such as the most common RANS models and large eddy
simulation with respect to drag predictions for cycling aerodynamics. They
found the $k$-$\omega$ SST model to perform best in comparison with
experimental data (i.e. force and moment areas) from wind tunnel measurements,
with average discrepancies of $6\%$ and consistently below $11\%$.

\begin{overpic}[width=390.25534pt]{Figure1.pdf} \end{overpic} Figure 1:
Geometry of a formation of runners surrounded by the computational domain and
the imposed boundary conditions (a). The axial spacing of $1.2\,m$ and the
lateral spacing of $0.7\,m$ between the runners are indicated (b).

Moreover, Crouch et al. [2016] has shown for cycling aerodynamics that the
effect of the dynamical motion of the limbs on instantaneous drag is only
minor. Additionally, wind tunnel experiments by Inoue et al. [2016] revealed
an approximately $10\%$ increase in drag values of solo running when using a
moving-belt system as compared with a stationary setup. In their measurements
the formation of unrealistic ground boundary layers was thus inhibited, as it
is the case in this study. We use a stationary model of an athlete together
with a moving ground boundary, which was found to be sufficient for the
prediction of the global aerodynamic forces.

Furthermore, discretization in this study is performed with a second-order
scheme for the spatial derivatives. Pressure and velocity are calculated in a
coupled approach. The fluid is defined as air at ambient temperature of
$T=$293.15\text{\,}\mathrm{K}$=$20\text{\,}\mathrm{\SIUnitSymbolCelsius}$$ and
a kinematic viscosity of
$\nu=$1.516\text{\times}{10}^{-5}\text{\,}{\mathrm{m}}^{2}\text{\,}{\mathrm{s}}^{-1}$$.
Additionally, the velocity inlet of the domain pictured in Fig. 1(a) imposes
an axial flow velocity at the chosen running speeds, while the pressure outlet
applies an ambient pressure of
$p=$101\,325\text{\,}\mathrm{Pa}$=$1\text{\,}\mathrm{atm}$$. The ground is
defined as wall moving at the running speed, which ensures that no boundary
layers form at the bottom that could invalidate the results. A Dirichlet
boundary condition, $u_{i}=0$, is applied for all components of the flow
velocity directly at the ground. The surrounding walls are symmetry boundaries
with frictionless properties. This is achieved through a Neumann boundary
condition $\partial u_{i}/\partial x_{i}=0$.

Tab. 1 gives the various speeds of this study with the respective Reynolds
numbers of the flow, using the runner's height for the definition of the
characteristic length. While all four speeds are considered for the assessment
of single-runner aerodynamics (cf. Sec. 3.1), the drafting formations are
studied at a running speed of
$21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ (cf. Sec. 3.2-3.4). This is
a high speed occurring during elite marathon running and allows for good
comparison between the considered formations. In fact, the current official
marathon world record for men is 2:01:39 by Eliud Kipchoge (September 16, 2018
during the Berlin Marathon; at an average speed of
$\sim$21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$$) and for women is
2:14:04 by Brigid Kosgei (October 13, 2019 during the Chicago Marathon; at an
average speed of $\sim$19\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$$).
Here, we chose a running velocity at the upper range of these numbers to
assess the ideal potential of drafting with respect to metabolic quantities,
drafted velocity, and time savings.

Table 1: Running speeds and associated Reynolds numbers $Re=\left(U_{\infty}L\right)/\nu$ at a characteristic length of $L=$1.65\text{\,}\mathrm{m}$$. Running speed ($\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$) | Reynolds number
---|---
$15$ ($4.16\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) | $452\,770\text{\,}$
$18$ ($5\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) | $544\,195\text{\,}$
$21$ ($5.83\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) | $634\,531\text{\,}$
$36$ ($10\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) | $1\,088\,390\text{\,}$

The inherent challenge of simulations lies in ensuring that physical effects
are computed as accurately as possible. Here, that means that flow separation
lines and points, as well as boundary layer thickness and characteristics are
captured for the investigated Reynolds numbers, which lie in the turbulent
regime (cf. Tab. 1). In order to avoid inaccurate prediction of aerodynamic
forces on athletes at high Reynolds numbers, as described by Meile et al.
[2006], it is crucial to optimize near-surface grid resolutions, such that no
numerically induced jumps of separation points occur. This is important, as an
artificially delayed flow separation would lead to a sharp drop in the
pressure drag, which is the dominant component of total drag of bluff bodies
(such as the human body). A thorough verification and validation of the
applied numerical approach is presented in A. The solution of the discretized
continuity and momentum equations is performed with the solver _ANSYS Fluent
19.2_ , post-processing of flow variables is achieved using the post-
processing tool _ParaView_ and drag distribution and resulting power values
are calculated with the free scientific programming language _GNU Octave_.

###  2.2 Geometry and mesh

For the numerical setup we apply the realistic, three-dimensional geometry of
a female runner with a height of $L=$1.65\text{\,}\mathrm{m}$$, a projected
frontal area of $A\approx$0.44\text{\,}{\mathrm{m}}^{2}$$, and a body mass of
$m\approx$55\text{\,}\mathrm{kg}$$. The geometry represents an average of 59
individual track athletes, scanned in-house and positioned into a running
pose. The obtained pose-invariant statistical model for the female human body
was reconstructed from the high-quality scan data using the method outlined by
Colaianni et al. [2014]. The axial and lateral spacing between the runners is
chosen according to estimates of realistic running conditions. Polidori et al.
[2020] for instance uses an axial distance of $1.3\text{\,}\mathrm{m}$ and a
lateral shoulder-to-shoulder distance of $0.3\text{\,}\mathrm{m}$, which are
taken from an actual drafting situation during a marathon race.

The computational domain shown in Fig. 1(a) has a distance of
$7\text{\,}\mathrm{m}$ from the main runner to the upstream boundary and
$13\text{\,}\mathrm{m}$ to the downstream boundary. It contains far-field
refinement as well as three refinement stages close to the runner's surfaces.
Furthermore, 10 prism layers (i.e. inflation layers) are inserted directly at
all solid boundaries with a steady growth factor of the prism layer thickness
of $1.2$. These measures ensure that the dimensionless wall distance of
$y^{+}\approx 1$ holds in the entire domain and all boundary layers are well
resolved. Due to the local mesh refinement the cell count varies sharply
between formations depending on the number of runners included. It rises from
$40$ million cells for the single-runner case up to $107$ million cells for
the final formation including a total of four runners. The near-surface mean
cell size of $3.125\text{\times}{10}^{-3}\text{\,}\mathrm{m}$ and the outer-
boundary mean cell size of $0.196\text{\,}\mathrm{m}$ for the final grid
settings have been chosen following a grid convergence study (cf. A).
Homogeneous, gradual cell growth is applied between the near and far field.
Locally, especially in the boundary layer regions, perpendicular mesh
dimensions are significantly lower (more than by a factor 10) than the near-
surface mean cell size, which is the cell size of the finest grid level just
outside the prism layers. A depiction of the considered formations of this
study is shown in Fig. 4-5.

###  2.3 Calculation of power values

Various empirical models exist for the calculation of the mechanical power
output while running. Fukunaga et al. [1980] found through experiments on
athletic runners on force platforms a relation between running velocity $u$
and the exerted power for forward motion as

$P_{R,0}=0.436u^{2.01},$ (1)

with $u$ in $\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$ and $P_{R}$ in
$\text{\,}\mathrm{W}\text{\,}{\mathrm{kg}}^{-1}$. In a study using motion
tracking, Cavagna and Kaneko [1977] also took into account additional energy
expenditure created by the movement of the limbs. They arrive at a function
for the specific running power per unit mass as

$P_{R,0}=9.42+4.73u+0.266u^{1.993}.$ (2)

Here, $u$ is given in $\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ and
$P_{R,0}$ in
$\text{\,}\mathrm{cal}\text{\,}{\mathrm{kg}}^{-1}\text{\,}{\mathrm{min}}^{-1}$.
For our considered cases, both approaches give similar results with deviations
in the range of only $2$-$5\%$, however, Eq. (2) is used in the following
power calculations due to being the more complete model of running motion. By
`more complete' we mean that it takes into account both internal and external
work performed for the modelling of the total mechanical work

$W_{tot}=W_{ext}+W_{int}.$ (3)

Here, external work $W_{ext}$ is due to acceleration and lift of the centre of
mass, while internal work refers to both the translational and rotational
acceleration of the limbs relative to the trunk. This approach has also been
suggested more recently by the work of Saibene and Minetti [2003], Pavei et
al. [2019], and Gray et al. [2020].

Using Eq. (2), we arrive at the power

$P_{R}=P_{R,0}\cdot\frac{4.1868}{60}\cdot m$ (4)

in $W$ for a runner of mass $m$ in $\text{\,}\mathrm{kg}$, in the following
simply called _running power_. Furthermore, there is a certain aerodynamic
power $P_{A}$, used to overcome drag $F_{D}$, which is defined as

$P_{A}=F_{D}\cdot u=\frac{1}{2}\rho u^{3}C_{D}A.$ (5)

Here, the expression for $F_{D}$ assumes stagnant air and the absence of
crosswinds. As a result, the total mechanical power generated while running is
the sum of the running power $P_{R}$ and the aerodynamic power $P_{A}$:

$P_{tot}=P_{R}+P_{A}.$ (6)

Using the total mechanical power $P_{tot}$ and the associated efficiency
$\epsilon$ at a certain run speed by Cavagna and Kaneko [1977], the metabolic
power

$P_{meta}=\frac{P_{tot}}{\epsilon},$ (7)

which gives the metabolic energy expenditure during running, can be computed.

It must be emphasized here that the science behind efficiency of locomotion
and in particular of walking and running is far from concluded. In both
running and walking the muscles utilize energy stored during a previous phase
of stretching in the following phase of contraction. As stated by Cavagna and
Kaneko [1977], the efficiency of walking reaches a maximum of $0.35$-$0.40$ at
intermediate speeds due to the expected properties of the contractile
component of muscle. This maximum of efficiency follows the force-velocity
relation and the trends due to initial efficiency of muscle [Hill, 1964]. Also
in cycling the efficiency of locomotion was found to reach maximum values at
intermediate speeds ($\sim 0.22$ according to the early study of Dickinson
[1929], $\sim 0.25$ according to Ettema and Lorås [2009]). In running,
however, the efficiency increases steadily with speed (from $0.45$ to $0.70$),
suggesting that the positive work during that activity derives mainly from the
passive recoil of elastic muscle tissue and not from the active shortening of
contractile muscle components [Cavagna and Kaneko, 1977, Saunders et al.,
2004]. In support of this hypothesis, Cavagna et al. [1968] found that the
useful effect of pre-streched muscle on the performed work increases with the
speed of stretching and shortening. Furthermore, Komi and Bosco [1978]
investigated this utilization of elastic energy during eccentric-concentric
contraction of muscles in human motion and found that up to $90\%$ of energy
produced in the pre-stretching phase is recovered in case of counter-movement
jumps. More recently, Komi [2000] used data based on in-vivo force
measurements, buckle-transducer technique, and optic-fiber technique to show
that the stretch-shortening cycle in human skeletal muscle leads to
considerable performance enhancement compared to simple concentric action.
With respect to running economy, Hunter et al. [2015] showed a clear positive
correlation of stretch-shortening cycle potentiation and increases of running
economy and recommended eccentric force development (e.g. by resistance
training) for performance improvements. Nonetheless, we point out that
efficiency and total mechanical power models are inherently difficult to
define due to the sometimes conflicting understanding of internal and external
work performed during human locomotion and the question, whether internal work
should appear in total mechanical power and efficiency expressions, as has
been suggested by Winter [1979] and more recently by Minetti et al. [2001]
(see e.g. Ettema and Lorås [2009] for a discussion of these issues with
regards to cycling, and Williams [2000] with regards to running).

##  3 Results

###  3.1 Single-runner aerodynamics

Using the geometry of the single runner, the drag forces acting on various
body parts are computed over the range of relevant velocities (cf. Fig. 2(b)).
It can be seen that the drag acting on torso and legs is highest, with both
parts roughly affected by the same aerodynamic loads. Each make for a fraction
of $39$-$40\%$ of the total drag on the runner. Furthermore, arms are affected
considerably less, with a share of approximately $13\%$, and the head with
$8\%$ of the total drag. These estimates give an overview of the sensitivity
of the various parts of the body towards drag loads and are directly related
to the frontal area of each part (cf. Eq. (5)). Fig. 2(a) gives the total drag
on the runner for the four simulated running speeds. While the drag for lower
speeds of $15$-$21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$, which are
relevant for long distance disciplines like marathon, lies in the range of
$3.39$-$6.52\text{\,}\mathrm{N}$, it sharply rises to
$18.45\text{\,}\mathrm{N}$ at the higher considered speed of
$36\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$. This demonstrates the
(quadratically) increasing relevance of drag reduction at higher speeds.
Moreover, the associated power values are calculated using the definitions of
Sec. 2.3. Fig. 2(c) shows the share of aerodynamic and running power in the
total mechanical power over the considered range of velocity. However, the
actual energy per unit time, which needs to be made available during motion,
is the metabolic power (cf. Eq. (7)).

\begin{overpic}[width=216.81pt]{Figure2.pdf} \end{overpic} Figure 2: Total
drag over the range of considered running speeds (a), as well as drag acting
on the body parts of the runner (b). Total mechanical power, consisting of
running power computed according to the empirical model by Cavagna and Kaneko
[1977] and aerodynamic power, as a function of speed (c). Polynomial functions
of second order (aerodynamic forces and running power) and of third order
(aerodynamic and total mechanical power) are used for fitting the data.

Here, the aerodynamic component provides a fraction of the metabolic power of
$2.6\%$, $3.4\%$, $4.3\%$, and $8.5\%$ for
$4.16\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$
($15\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$),
$5\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$
($18\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$),
$5.83\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$
($21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$), and
$10\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$
($36\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$). These can be considered
as threshold values for the maximum savings in metabolic energy consumption
during drafted running.

The resulting wake regions are depicted in Fig. 3. They are visualized by
isosurfaces of total pressure, giving the boundary between the negative
pressure region in the wake and the positive pressure region due to recovery
of the flow. There is a notable decrease of the wake length from
$1.14\text{\,}\mathrm{m}$ at
$15\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ to
$0.98\text{\,}\mathrm{m}$ at
$36\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ due to increased turbulent
mixing. The wake length is a good indication of the most beneficial range for
drafting, as the suction provided by the negative total pressure in the wake
slightly pulls a following runner. However, even at a length larger than that
there will be a considerable drag-reducing effect due to flow deceleration of
the head wind.

\begin{overpic}[width=433.62pt]{Figure3.pdf} \end{overpic} Figure 3: Wake
region behind the runner visualized by isosurfaces of total pressure at
$p_{tot}=0$.

###  3.2 Drag comparison of formations

At first we consider the performance of the investigated formations with
respect to their effect on the main runner's drag at a running speed of
$21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$. It can be concluded from
Fig. 4 that formation $2$ leads to the largest reduction of drag of $75.6\%$,
followed by formation $1$ with $70.1\%$. Formations $3$ and $4$ cause
significantly less drag reduction of $41.3\%$ and $33.4\%$, respectively.

\begin{overpic}[width=216.81pt]{Figure4.pdf} \end{overpic} Figure 4: Sketch of
the investigated formations (a) together with their effect on the main
runner's drag (b). The formations are compared to the case of the single
runner. Axial and lateral spacing between the runners follows the values
introduced in Fig. 1(b).

Furthermore, considering the drag distribution across the whole body of the
main runner, it can be recognized that the dominant part of the drag acts on
legs and torso (cf. Sec. 3.1). Thus, it is crucial for a drafting formation to
reduce the drag in these areas. Looking at Fig. 5, it can be seen that
formations $1$-$2$ achieve this, while formations $3$-$4$ cannot adequately
shield the body of the main runner from oncoming air. Additionally, it is
shown that while the drag reduction on head and arms is fairly equal for
formations $1$-$3$, formation $4$ shows almost no improvement for the drag on
the arms. The reasons for the varying aerodynamic performances of formations
$1$-$4$ are elaborated further in Sec. 3.3.

\begin{overpic}[width=433.62pt]{Figure5.pdf} \end{overpic} Figure 5: Drag
acting on the various body parts of the main runner alongside the total drag
for all involved runners, including pacers, for the formations 1 (a), 2 (b), 3
(c), and 4 (d).

Using the results for formation $1$, the distance between main runner and
pacer was further varied to give the drag values shown in Fig. 6. An
exponential function of the type $y=a+b\cdot\exp\left(-c/x\right)$ with $a=0$,
$b=5.63$, and $c=1.28$ was used to fit the data and has the desired property
of asymptotically approaching zero drag at small distances and the single
runner's maximum drag at large distances. The values of the parameters were
estimated by least-squares optimization using the Levenberg-Marquardt
algorithm [Levenberg, 1944, Marquardt, 1963]. The absolute axial limit of
separation between the two runners, given the dimension and stride length of
the applied geometry, was found to be approximately $0.6\text{\,}\mathrm{m}$.
In reality, however, this limit will barely be reached and the runners will be
unable to approach already considerably earlier. Pugh [1970] estimates a
natural limit of axial spacing of $1\text{\,}\mathrm{m}$ for middle-distance
running and states a possible drag reduction at this distance of up to $80\%$.
In this study we find a drag of $1.45\text{\,}\mathrm{N}$ acting on the main
runner at a distance of $1\text{\,}\mathrm{m}$, which corresponds to $22.2\%$
of the drag on the single runner and thereby lies in good agreement with the
results by Pugh [1970].

\begin{overpic}[width=216.81pt]{Figure6.pdf} \end{overpic} Figure 6: Drag on
the main runner as a function of the distance to the leading pacer in
formation 1. The numerical data is fitted using an exponential function with a
determination coefficient of $R^{2}=0.99$. In reality, there is an axial limit
which prohibits the further approach of the main runner to the pacer. This
limit is estimated from geometric considerations as $0.6\text{\,}\mathrm{m}$.
Pugh [1970] estimates this minimum distance during drafting at
$1.0\text{\,}\mathrm{m}$.

###  3.3 Aerodynamic effects on runners within formations

Air resistance while running manifests itself in static pressure increase as
the air is being decelerated and diverted at the front of the runner. Fig. 7
provides insight to the magnitude of the pressure acting on the main runner
for the various formations. The pressure coefficient

$c_{p}=\frac{p-p_{\infty}}{\frac{1}{2}\rho
U_{\infty}^{2}}=\frac{p-p_{\infty}}{p_{0}-p_{\infty}},$ (8)

with the free-stream pressure of
$p_{\infty}=$101\,325\text{\,}\mathrm{Pa}$=$1\text{\,}\mathrm{atm}$$ and a
relative free-stream velocity of
$U_{\infty}=$5.83\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$=$21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$$,
indicates regions of high or low relative pressure.

\begin{overpic}[width=433.62pt]{Figure7.pdf} \end{overpic} Figure 7: Pressure
coefficient on the main runner for the single-runner case (a) compared to the
formations 1-4 (b)-(e).

Fig. 7 confirms the beneficial pressure distribution for formations $1$ and
$2$. This means that there are no large areas of high pressure coefficient at
the front of the runner, nor any large regions of low negative pressure
coefficient at the back, which would be an indication for suction force
slightly pulling the runner back. Formation $3$ on the other hand shows a high
positive pressure coefficient at the front, which already hints at the leakage
of oncoming flow between the front pacers that can be seen in Fig. 8(c) and
Fig. 9(c). While formation $4$ does not display any significantly high
pressure peaks at the front, it is the elevated negative pressure coefficient
at the back that makes it less advantageous than formations $1$-$2$ (cf. Sec.
3.2).

\begin{overpic}[width=433.62pt]{Figure8.pdf} \end{overpic} Figure 8:
Streamlines seeded at the front of the leading pacers showing the diversion of
air around the main runner for all four formations. Streamlines are coloured
by velocity magnitude and the runners' surfaces by pressure coefficient.
\begin{overpic}[width=433.62pt]{Figure9.pdf} \end{overpic} Figure 9: Velocity
magnitude contours in the horizontal plane section through the runners'
centre. Formations 1-4 are shown in images (a)-(d).

Fig. 8-9 both demonstrate the reasons for the large differences in aerodynamic
performance between the formations. In case of formations $1$ and $2$ the main
runner is well embedded inside the wake region of the front pacer and
therefore only affected by low-speed, almost stagnant air relative to herself.
Formation $2$ has the additional advantage of a rear pacer, whose high-
pressure stagnation region at the front positively affects the low-pressure
wake region of the main runner by raising the pressure coefficient to a value
closer to zero, thus lowering suction. This can be seen when comparing Fig.
10(a) and Fig. 10(b). Furthermore, Fig. 8(c) and Fig. 9(c) clearly demonstrate
a major weakness of formation $3$, which is the passing air between the front
pacers. This leads to a large region of high stagnation pressure directly at
the area of the impacting air on the main runner, which is also visible in
Fig. 7(d) and Fig. 10(c). Also formation $4$ is suboptimal due to air being
accelerated between the main runner and the pacers on the side. This leads to
the aforementioned effect of higher drag at the arms compared to the other
formations (cf. Sec. 3.2). Additionally, this effect of locally increased flow
velocity between the runners leads to a shortening of the wake region of slow
air for the main runner and thus negatively impacts her drafting conditions
(cf. Fig. 9(d)).

\begin{overpic}[width=325.215pt]{Figure10.pdf} \end{overpic} Figure 10:
Pressure coefficient in the vertical plane section through the main runner's
centre. Formations 1-4 are depicted in images (a)-(d).

###  3.4 Metabolic power savings and performance predictions

Using the drag values from Sec. 3.2, we can compute the associated aerodynamic
power and the total mechanical power using Eq. (2)-(6). Tab. 2 lists the
results for the single runner and the four formations. Additionally, relative
changes in aerodynamic power and total power are given by comparison with the
results of the single-runner case. The resulting metabolic power needed for
running in the considered formations is computed by using the efficiency of
running at the given velocity as detailed in Cavagna and Kaneko [1977].
Furthermore, running economy as the specific metabolic rate per unit mass is
given together with its relative changes in Tab. 3. As described by Kipp et
al. [2019], the relationship between improvements in running economy and
improvements in velocity is not a linear, but rather a curvilinear one. This
means that for a given improvement in running economy, velocity gains will be
smaller at higher speeds and larger at lower speeds. At elite marathon speeds
of around $5.5\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$ there is roughly
a $2/3^{\text{rds}}$ increase in velocity for each percent increase in running
economy. By applying the conversion factors of Kipp et al. [2019] we arrive at
estimates for the velocity gains for each formation in Tab. 4. The relative
changes in speed are $2.1\%$ (formation 1), $2.3\%$ (formation 2), $1.3\%$
(formation 3), and $1.1\%$ (formation 4).

Table 2: Aerodynamic power $P_{A}$, and total mechanical power $P_{tot}$ computed with Eq. (5) and (6) and by computing the running power $P_{R}=$858.1\text{\,}\mathrm{W}$$ at $21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ ($5.83\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) using the model of Cavagna and Kaneko [1977]. Savings $\Delta P_{A}$ and $\Delta P_{tot}$ are given with respect to the case of the runner alone. Formation |  $F_{D}$ ($\text{\,}\mathrm{N}$) |  $P_{A}$ ($\text{\,}\mathrm{W}$) | $\Delta P_{A}$ |  $P_{tot}$ ($\text{\,}\mathrm{W}$) | $\Delta P_{tot}$
---|---|---|---|---|---
Single runner | $6.52$ | $38.0$ | - | $896.1$ | -
$1$ | $1.96$ | $11.4$ | $70.0\%$ | $869.5$ | $3.0\%$
$2$ | $1.58$ | $9.2$ | $75.8\%$ | $867.3$ | $3.2\%$
$3$ | $3.84$ | $22.4$ | $41.1\%$ | $880.5$ | $1.7\%$
$4$ | $4.33$ | $25.3$ | $33.4\%$ | $883.4$ | $1.4\%$
Table 3: Metabolic power $P_{meta}$ according to Eq. (7) and running economy $RE$, as well as percent improvements of the running economy $\Delta RE$ for each of the formations with respect to the single-runner case. Efficiency of locomotion for a running velocity of $21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ ($5.83\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) is taken as $\epsilon=0.64$ [Cavagna and Kaneko, 1977]. Formation |  $P_{meta}$ ($\text{\,}\mathrm{W}$) |  $RE$ ($\text{\,}\mathrm{W}\text{\,}{\mathrm{kg}}^{-1}$) | $\Delta RE$
---|---|---|---
Single runner | $1400.2$ | $25.5$ | -
1 | $1358.6$ | $24.7$ | $3.1\%$
2 | $1355.2$ | $24.6$ | $3.5\%$
3 | $1375.8$ | $25.0$ | $2.0\%$
4 | $1380.3$ | $25.1$ | $1.6\%$
Table 4: Drafted running speeds $u_{draft}$ for all formations alongside potential time $t_{draft}$ and time savings through drafting $\Delta t_{draft}$ over the marathon distance of $42\,195\text{\,}\mathrm{m}$ by comparison with the undrafted single-runner case. Formation |  $u_{draft}$ ($\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) |  $t_{draft}$ ($\text{\,}\mathrm{s}$) |  $\Delta t_{draft}$ ($\text{\,}\mathrm{s}$)
---|---|---|---
1 | $5.95$ |  $7092$ ($\approx$1.970\text{\,}\mathrm{h}$$) |  $142$ ($\approx$2.4\text{\,}\mathrm{min}$$)
2 | $5.96$ |  $7080$ ($\approx$1.967\text{\,}\mathrm{h}$$) |  $154$ ($\approx$2.6\text{\,}\mathrm{min}$$)
3 | $5.91$ |  $7140$ ($\approx$1.983\text{\,}\mathrm{h}$$) |  $98$ ($\approx$1.6\text{\,}\mathrm{min}$$)
4 | $5.89$ |  $7164$ ($\approx$1.990\text{\,}\mathrm{h}$$) |  $70$ ($\approx$1.2\text{\,}\mathrm{min}$$)

##  4 Discussion and conclusions

Drafting and the resulting aerodynamic drag reduction cause a notable effect
on running and the energy expenditure used for it. This study shows the
overall aerodynamics and associated energy expenditure during running for a
range of relevant speeds. Considerable shares of up to $8.5\%$ at a velocity
of $36\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ are computed.

Additionally, differences in drag, pressure distribution across the runners,
and potential power savings are demonstrated for four formations that could
arise during running competitions. As shown in Sec. 3.2, simple formations of
axial spacing of the pacers (i.e. one pacer in front, as in formation 1, or
one pacer in front and one in the back, as in formation 2) already lead to
considerable reduction of over $70\%$ in drag at a distance of
$1.2\text{\,}\mathrm{m}$. Furthermore, formations 1 and 2 also show a clear
positive effect on leading and trailing pacers, since their overpressure and
underpressure regions interact with the main runner. This effect has also been
noted by Beaumont et al. [2019] together with a measurable advantage for the
pacers with respect to oxygen consumption.

We further deduce an empirical relationship of drag force acting on the main
runner versus her distance to a leading pacer, which can be used for future
prediction of a pacer's aerodynamic effects (cf. Fig. 6). The results are
congruent with experimental data by Pugh [1970], who predicts drag decrease of
up to $80\%$ when running behind a pacer at $1\text{\,}\mathrm{m}$. Davies
[1980], who claimed possible drag savings in the range of $80$-$85\%$, also
confirms these results.

Moreover, the drag distribution at various body parts is examined in Sec. 3.2.
From this and the explanations in Sec. 3.3, we can see a large sensitivity of
the drag acting on the main runner if pacers are positioned laterally, such as
in formations $3$ and $4$. In formation $3$, flow leaks between the front
pacers, causing large pressure peaks at the runner's centre (cf. Fig. 7).
Formation $4$ on the other hand has the negative effect of pacers at the
runner's sides and the resulting acceleration of air between the runners leads
to a relatively low drag reduction of only $33.4\%$ compared to formation $1$
with $70.1\%$. This suggests that runners at the side should be avoided if
optimal drag savings are the goal.

In essence, there are three aerodynamic effects of drafting that lead to a
reduction of axial force (i.e. drag) while running:

  1. (i)

_Shielding against fast, oncoming air._
Therefore, there is no sharp deceleration of air and resulting formation of
areas of high stagnation pressure on the runner's surface.

  2. (ii)

_Suction effects due to the region of negative total pressure in the leading
runner's wake._
Any object moving through air forms a wake region of relatively low pressure.
This underpressure can be harnessed to the benefit of a trailing runner. The
wake region of negative total pressure has been found in this study to extend
to about $1.14\text{\,}\mathrm{m}$ at
$15\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ down to
$0.98\text{\,}\mathrm{m}$ at
$36\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ from the trailing runner's
back.

  3. (iii)

_Delay of separation points due to turbulent free shear layers._
Turbulence levels are naturally increased in the high-shear-stress regions of
the wake, where air is displaced against the stagnant surroundings. This
raises momentum transfer in the boundary layers of the runner and delays
separation due to formation of a turbulent boundary layer from a laminar one.
As a result, pressure drag (i.e. form drag) drops, which is the dominant
source of drag of bluff bodies.

Eventually, the aerodynamic power savings are computed in this study and
listed together with the total mechanical power and the metabolic power output
in Sec. 3.4. The derived metabolic power and running economy of up to
$1400.2\text{\,}\mathrm{W}$ and
$25.5\text{\,}\mathrm{W}\text{\,}{\mathrm{kg}}^{-1}$ for the single runner at
a speed of $21\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$ agrees well with
values from literature for trained athletes [Kipp, 2017, Batliner et al.,
2018]. The largest improvements of running economy occur for formation $1$
with $3.1\%$ and formation $2$ with $3.5\%$. These metabolic gains are
slightly higher than the $2.84\%$ gain computed by Polidori et al. [2020] for
the best position of Kenenisa Bekele's drafting strategy during the 2019
Berlin marathon. The reason for this is most likely the larger axial spacing
of $1.3\text{\,}\mathrm{m}$ to the front pacer as opposed to
$1.2\text{\,}\mathrm{m}$ in the current study. Furthermore, a particular
benefit of this study is the additional performance prediction in terms of
improvements of velocity and time for the investigated formations. By applying
the curvilinear relationship between gains in velocity and running economy
suggested by Kipp et al. [2019], we deduce possible speed gains in the range
of $1.1$-$2.3\%$. Time savings of $1.2\text{\,}\mathrm{min}$ (formation 4) to
$2.6\text{\,}\mathrm{min}$ (formation 2) are achieved over the marathon
distance.

In this study we focus on the impact of aerodynamic drafting on running
economy. However, there are a multitude of other factors that influence
running economy and performance. With respect to biomechanics, some of these
factors are (lower) limb mass distribution, Achilles tendon moment arm and
lower body musculotendinous structures for reutilization of elastic energy, as
well as running style and gait patterns (e.g. so-called _pose running_ to
benefit from elastic recoil forces in the lower limbs), among others [Williams
and Cavanagh, 1987, Saunders et al., 2004, Barnes and Kilding, 2015].

Concluding from this work, it is possible to give targeted suggestions for
formation patterns that are best suited for cooperative drafting. It is shown
that axially positioned pacers give the best performance and that laterally
positioned pacers should be running at a distance further away from the main
runner to avoid local flow acceleration, while still contributing positively
to drag reduction. Such an approach has been applied with the recent
successful world record attempt by Eliud Kipchoge during the _INEOS 1:59
Challenge_ , where a staggered formation in the shape of an inverted V has
been chosen. Future studies could explore a larger number of formations by CFD
and subsequently derive reduced-order models for parameters such as number of
pacers, or axial and lateral spacing.

Acknowledgements
This work was realised through research budget from ADIDAS and KTM
E-TECHNOLOGIES. Computational resources were provided by the parallel
computing cluster at KTM in Mattighofen, Austria.

Conflict of interest statement
The authors have nothing to disclose that would have biased this work.

##  Appendix A Verification and validation

Since errors and uncertainties are an unavoidable part of any simulation
method such as CFD, measures must be taken to quantify their levels and to
minimize them. According to the guidelines by the American Institute of
Aeronautics and Astronautics [1998], both verification and validation form a
crucial first step of any CFD simulation. While verification can be seen as
the process of ensuring that `the equations are solved right', validation
means ensuring `that the right equations are solved' [Roache, 1998].
In this work, we first perform a grid convergence study by successively
refining the mesh until no significant changes in the relevant quantities are
occurring at the chosen high-order schemes (i.e. verification). Then we
compare the obtained drag coefficients of typical bluff bodies to well-
established values from literature. Finally, since no wind tunnel measurement
results are available for the specific runner geometry we use, we compare the
numerically computed drag coefficient of the runner to experimental values of
previous studies (i.e. validation).

#### Grid convergence study

The average cell size of the mesh at the finest grid level outside the prism
layers has been refined by a factor of two to obtain the target cell sizes at
the runner's surface as illustrated in Fig. 11. All grids contain $5$
inflation layers (i.e. prism layers) at the object surfaces, except the final
one, which has $10$ inflation layers. It can be seen in Fig. 11 that for both
the drag of the considered bluff bodies as well as for the applied runner
geometry the resulting changes in the final refinement stages are small.
Following Roache [1997], we can calculate the discretization error of two
meshes with refinement ratio $r=h_{2}/h_{1}=0.5$ and with cell sizes $h_{1}$
and $h_{2}$ by using the difference $\phi_{2}-\phi_{1}$ between the two
solutions of the target quantity, as

$\displaystyle E_{\phi,1}$ $\displaystyle=$
$\displaystyle\frac{\phi_{2}-\phi_{1}}{1-r^{p}},$ (9) $\displaystyle
E_{\phi,2}$ $\displaystyle=$ $\displaystyle
r^{p}\left(\frac{\phi_{2}-\phi_{1}}{1-r^{p}}\right),$ (10)

with the order of the numerical scheme $p=2$. Using this procedure, we obtain
$E_{\phi,1}=[E_{\phi,1,cube},E_{\phi,1,cylinder},E_{\phi,1,sphere}]=[0.0027,0.0213,0.0240]$
and
$E_{\phi,2}=[E_{\phi,2,cube},E_{\phi,2,cylinder},E_{\phi,2,sphere}]=[0.0007,0.0053,0.0060]$
for the cell sizes $0.006\,25\text{\,}\mathrm{m}$ and
$0.003\,12\text{\,}\mathrm{m}$, both with five inflation layers. For the
applied runner geometry the discretization errors are computed as
$E_{\phi,1,runner}=0.0027$ and $E_{\phi,2,runner}=0.0007$. Already with the
coarser mesh the error is consistently below $3\%$ for all considered bluff
bodies, and with the fine mesh it is below $1\%$. Apart from the
discretization error there is also a linearization error inherent to all CFD
simulations, which comes from the linearization of the governing equations
through numerical schemes. Minimization of the linearization errors of the
simulations is ensured through monitoring of the residuals of all flow
variables and the convergence of the aerodynamic forces. We apply a second-
order scheme for increased accuracy as well as tolerance values for all
residuals of $\mathcal{O}(10^{-6})$.

\begin{overpic}[width=216.81pt]{Figure11.pdf} \end{overpic} Figure 11: Mesh
convergence demonstrated by the final four refinement stages. The drag
coefficient of typical bluff bodies (a) and of the runner geometry used for
this study (b) are shown (i.l. stands for inflation layers).

#### Drag coefficients of general bluff bodies

Using the numerical results for the drag coefficients of the bluff bodies, a
comparison is done with experimental data from literature. The results are
summarized in Tab. 5.

Table 5: Numerical values from this study and experimental values from literature for typical bluff-body drag coefficients $C_{D}$ at a considered free-stream flow velocity of $U_{\infty}=$18\text{\,}\mathrm{km}\text{\,}{\mathrm{h}}^{-1}$$ ($5\text{\,}\mathrm{m}\text{\,}{\mathrm{s}}^{-1}$) and a Reynolds number of $Re\approx(1.5\text{-}1.6)\times 10^{5}$. Object |  $C_{D}$ (num.) |  $C_{D}$ (exp.)
---|---|---
Cube | $1.06$ |  $1.05$ [Hoerner, 1951]
Cylinder | $0.58$ |  $0.60$ [Prosser and Smith, 2015]
Sphere | $0.38$ |  $0.40$ [Munson et al., 1990]

It can be seen that there is an overall good agreement between the numerical
results and literature, despite the fact that the considered Reynolds number
is close to the drag crisis for the sphere and the cylinder in crossflow. In
this range there is a transition from laminar to turbulent boundary layers,
which leads to a sudden drop in the viscous drag due to delayed flow
separation. It is therefore crucial to properly resolve the flow field across
the boundary layer, which is achieved here through insertion of a sufficient
number of prism layers.

#### Drag coefficients of runners

In the final part of the validation the runner's drag is compared to the
results of several experimental studies on runners.

Table 6: Comparison of values from literature in chronological order for the drag coefficient $C_{D}$ of humans at various activities. Authors | Activity | $C_{D}$
---|---|---
Walpert and Kyle [1989] | Running | $0.79$
Davies [1980] | Running |  $0.82$-$0.91$
Pugh [1971] | Running | $0.8$
| Walking | $0.7$
Hoerner [1965] | Standing |  $1.0$-$1.3$ (clothed)
|  |  $0.9$-$1.2$ (nude)
Hill [1928] | Running | $0.9$
| Standing | $0.98$

From the data in Tab. 6 it is obvious that a spread of the drag coefficient
values exists with respect to the different activities, but also for running
itself. It proves the difficulty of establishing well-defined drag
coefficients for humans experimentally. However, there is a clear reduction of
the drag coefficient from standing to walking and running. This has been
argued to be due to the limbs being partially lifted in the air leading to a
temporary reduction of frontal area and the subject being less affected by the
air resistance over time. Furthermore, clothing has a significant impact as
well, leading to increased drag and possibly earlier flow separation,
especially if loose [Nørstrud, 2009]. The quantitative effect of surface
properties on running aerodynamics has been assessed by Kyle and Caiozzo
[1986], who estimate the increase of drag by loose clothing to be up to
$4.2\%$ and by hair to be $4$-$6\%$, depending on length. In the current
study, we obtain a drag coefficient of the considered runner of $0.74$ at a
comparable Reynolds number, which is close to the values established by Pugh
[1971], Davies [1980], and Walpert and Kyle [1989]. One reason, why the drag
coefficient computed in our case is slightly lower than the values from
literature given in Tab. 6 might be the smooth geometry, allowing for delayed
flow separation and thus reduced form drag on the runner. This is an effect
that would not occur in a wind tunnel setup involving real runners or their
scaled models due to the aforementioned conditions of surface roughness, hair,
and clothing.

## References

  * American Institute of Aeronautics and Astronautics [1998] American Institute of Aeronautics and Astronautics .  AIAA Guide for the Verification and Validation of Computational Fluid Dynamics Simulations, 1998.
  * Barnes and Kilding [2015] Barnes KR, Kilding AE.  Running economy: Measurement, norms, and determining factors.  Sports Medicine 2015;1(1):1–15.
  * Bassett Jr et al. [1991] Bassett Jr DR, Flohr J, Duey WJ, Howley ET, Pein RL.  Metabolic responses to drafting during front crawl swimming.  Medicine and Science in Sports and Exercise 1991;23(6):744–7.
  * Batliner et al. [2018] Batliner ME, Kipp S, Grabowski AM, Kram R, Byrnes WC.  Does metabolic rate increase linearly with running speed in all distance runners?  Sports Medicine International Open 2018;2(1):E1.
  * Beaumont et al. [2019] Beaumont F, Bogard F, Murer S, Polidori G, Madaci F, Taiar R.  How does aerodynamics influence physiological responses in middle-distance running drafting?  Mathematical Modelling of Engineering Problems 2019;6:129–35.
  * Blocken et al. [2013] Blocken B, Defraeye T, Koninckx E, Carmeliet J, Hespel P.  CFD simulations of the aerodynamic drag of two drafting cyclists.  Computers & Fluids 2013;71:435–45.
  * Blocken et al. [2018] Blocken B, Toparlar Y, van Druenen T, Andrianne T.  Aerodynamic drag in cycling team time trials.  Journal of Wind Engineering and Industrial Aerodynamics 2018;182:128–45.
  * Cavagna et al. [1968] Cavagna G, Dusman B, Margaria R.  Positive work done by a previously stretched muscle.  Journal of Applied Physiology 1968;24(1):21–32.
  * Cavagna and Kaneko [1977] Cavagna G, Kaneko M.  Mechanical work and efficiency in level walking and running.  Journal of Physiology 1977;268(2):467–81.
  * Chatard and Wilson [2003] Chatard JC, Wilson B.  Drafting distance in swimming.  Medicine and Science in Sports and Exercise 2003;35(7):1176–81.
  * Colaianni et al. [2014] Colaianni M, Zollhöfer M, Süßmuth J, Seider B, Greiner G.  A pose invariant statistical shape model for human bodies.  In: Proceedings of the 5th International Conference on 3D Body Scanning Technologies. 2014. p. 327–36.
  * Crouch et al. [2016] Crouch TN, Burton D, Thompson MC, Brown NA, Sheridan J.  Dynamic leg-motion and its effect on the aerodynamic performance of cyclists.  Journal of Fluids and Structures 2016;65:121–37.
  * Davies [1980] Davies C.  Effects of wind assistance and resistance on the forward motion of a runner.  Journal of Applied Physiology 1980;48(4):702–9.
  * Defraeye et al. [2010] Defraeye T, Blocken B, Koninckx E, Hespel P, Carmeliet J.  Computational fluid dynamics analysis of cyclist aerodynamics: Performance of different turbulence-modelling and boundary-layer modelling approaches.  Journal of Biomechanics 2010;43(12):2281–7.
  * Dickinson [1929] Dickinson S.  The efficiency of bicycle-pedalling, as affected by speed and load.  Journal of Physiology 1929;67(3):242–55.
  * Edwards and Byrnes [2007] Edwards AG, Byrnes WC.  Aerodynamic characteristics as determinants of the drafting effect in cycling.  Medicine and Science in Sports and Exercise 2007;39(1):170–6.
  * Ettema and Lorås [2009] Ettema G, Lorås HW.  Efficiency in cycling: A review.  European Journal of Applied Physiology 2009;106(1):1–14.
  * Fitton et al. [2018] Fitton B, Caddy O, Symons D.  The impact of relative athlete characteristics on the drag reductions caused by drafting when cycling in a velodrome.  Proceedings of the Institution of Mechanical Engineers, Part P: Journal of Sports Engineering and Technology 2018;232(1):39–49.
  * Fukunaga et al. [1980] Fukunaga T, Matsuo A, Yuasa K, Fujimatsu H, Asahina K.  Effect of running velocity on external mechanical power output.  Ergonomics 1980;23(2):123–36.
  * Fuss [2018] Fuss FK.  Slipstreaming in gravity powered sports: Application to racing strategy in ski cross.  Frontiers in Physiology 2018;9:1032.
  * Gray et al. [2020] Gray A, Andrews M, Waldron M, Jenkins D.  A model for calculating the mechanical demands of overground running.  Sports Biomechanics 2020;:1–22.
  * Hellsten [1998] Hellsten A.  Some improvements in Menter's k-$\omega$ SST turbulence model.  In: 29th AIAA, Fluid Dynamics Conference. 1998\. p. 2554.
  * Hill [1928] Hill AV.  The air-resistance to a runner.  Proceedings of the Royal Society of London Series B, Containing Papers of a Biological Character 1928;102(718):380–5.
  * Hill [1964] Hill AV.  The efficiency of mechanical power development during muscular shortening and its relation to load.  Proceedings of the Royal Society of London Series B Biological Sciences 1964;159(975):319–24.
  * Hoerner [1951] Hoerner SF.  Aerodynamic Drag: Practical Data on Aerodynamic Drag Evaluated and Presented.  Otterbein Press, 1951.
  * Hoerner [1965] Hoerner SF.  Fluid-Dynamic Drag: Theoretical, experimental and statistical information.  SF Hoerner Fluid Dynamics, 1965.
  * Hoogkamer et al. [2018] Hoogkamer W, Snyder KL, Arellano CJ.  Modeling the benefits of cooperative drafting: Is there an optimal strategy to facilitate a sub-2-hour marathon performance?  Sports Medicine 2018;48(12):2859–67.
  * Hunter et al. [2015] Hunter GR, McCarthy JP, Carter SJ, Bamman MM, Gaddy ES, Fisher G, Katsoulis K, Plaisance EP, Newcomer BR.  Muscle fiber type, achilles tendon length, potentiation, and running economy.  Journal of Strength & Conditioning Research 2015;29(5):1302–9.
  * Inoue et al. [2016] Inoue T, Okayama T, Teraoka T, Maeno S, Hirata K.  Wind-tunnel experiment on aerodynamic characteristics of a runner using a moving-belt system.  Cogent Engineering 2016;3(1):1231389.
  * Katz [1995] Katz J.  Race Car Aerodynamics: Designing for Speed.  Bentley Publishers, 1995.
  * Kipp [2017] Kipp S.  Why does metabolic rate increase curvilinearly with running velocity?  Master's thesis; University of Colorado, Boulder, Colorado; 2017\.
  * Kipp et al. [2019] Kipp S, Kram R, Hoogkamer W.  Extrapolating metabolic savings in running: Implications for performance predictions.  Frontiers in Physiology 2019;10:79.
  * Komi [2000] Komi PV.  Stretch-shortening cycle: A powerful model to study normal and fatigued muscle.  Journal of Biomechanics 2000;33(10):1197–206.
  * Komi and Bosco [1978] Komi PV, Bosco C.  Utilization of stored elastic energy in leg extensor muscles by men and women.  Medicine and Science in Sports 1978;10(4):261–5.
  * Kyle and Caiozzo [1986] Kyle CR, Caiozzo VJ.  The effect of athletic clothing aerodynamics upon running speed.  Medicine and Science in Sports and Exercise 1986;18(5):509–15.
  * Levenberg [1944] Levenberg K.  A method for the solution of certain non-linear problems in least squares.  Quarterly of Applied Mathematics 1944;2(2):164–8.
  * Malizia and Blocken [2020] Malizia F, Blocken B.  Bicycle aerodynamics: History, state-of-the-art and future perspectives.  Journal of Wind Engineering and Industrial Aerodynamics 2020;200:104134.
  * Marquardt [1963] Marquardt DW.  An algorithm for least-squares estimation of nonlinear parameters.  Journal of the Society for Industrial and Applied Mathematics 1963;11(2):431–41.
  * Meile et al. [2006] Meile W, Reisenberger E, Mayer M, Schmölzer B, Müller W, Brenn G.  Aerodynamics of ski jumping: Experiments and CFD simulations.  Experiments in Fluids 2006;41(6):949–64.
  * Menter [1992] Menter FR.  Improved two-equation k-$\omega$-turbulence models for aerodynamic flows.  National Aeronautics and Space Administration, Ames Research Center, Moffett Field, CA 1992;.
  * Menter [1994] Menter FR.  Two-equation eddy-viscosity turbulence models for engineering applications.  AIAA Journal 1994;32(8):1598–605.
  * Minetti et al. [2001] Minetti AE, Pinkerton J, Zamparo P.  From bipedalism to bicyclism: Evolution in energetics and biomechanics of historic bicycles.  Proceedings of the Royal Society of London Series B: Biological Sciences 2001;268(1474):1351–60.
  * Munson et al. [1990] Munson BR, Young DF, Okiishi TH.  Fundamentals of Fluid Mechanics.  John Wiley & Sons, 1990.
  * Nørstrud [2009] Nørstrud H.  Sport Aerodynamics.  Springer, 2009.
  * Pavei et al. [2019] Pavei G, Zamparo P, Fujii N, Otsu T, Numazu N, Minetti AE, Monte A.  Comprehensive mechanical power analysis in sprint running acceleration.  Scandinavian Journal of Medicine & Science in Sports 2019;29(12):1892–900.
  * Polidori et al. [2020] Polidori G, Legrand F, Bogard F, Madaci F, Beaumont F.  Numerical investigation of the impact of Kenenisa Bekele's cooperative drafting strategy on its running power during the 2019 Berlin marathon.  Journal of Biomechanics 2020;107:109854.
  * Prosser and Smith [2015] Prosser D, Smith M.  Aerodynamics of finite cylinders in quasi-steady flow.  In: 53rd AIAA Aerospace Sciences Meeting. 2015\. p. 1931.
  * Pugh [1970] Pugh LGE.  Oxygen intake in track and treadmill running with observations on the effect of air resistance.  Journal of Physiology 1970;207(3):823–35.
  * Pugh [1971] Pugh LGE.  The influence of wind resistance in running and walking and the mechanical efficiency of work against horizontal or vertical forces.  Journal of Physiology 1971;213(2):255–76.
  * Roache [1997] Roache PJ.  Quantification of uncertainty in computational fluid dynamics.  Annual Review of Fluid Mechanics 1997;29(1):123–60.
  * Roache [1998] Roache PJ.  Verification and Validation in Computational Science and Engineering.  Hermosa Publishers, 1998.
  * Romberg et al. [1971] Romberg G, Chianese F, Lajoie R.  Aerodynamics of race cars in drafting and passing situations.  Technical Report; SAE; 1971.
  * Rundell [1996] Rundell KW.  Effects of drafting during short-track speed skating.  Medicine and Science in Sports and Exercise 1996;28(6):765–71.
  * Saibene and Minetti [2003] Saibene F, Minetti AE.  Biomechanical and physiological aspects of legged locomotion in humans.  European Journal of Applied Physiology 2003;88(4-5):297–316.
  * Saunders et al. [2004] Saunders PU, Pyne DB, Telford RD, Hawley JA.  Factors affecting running economy in trained distance runners.  Sports Medicine 2004;34(7):465–85.
  * Walpert and Kyle [1989] Walpert RA, Kyle CR.  Aerodynamics of the human body in sports.  Journal of Biomechanics 1989;22(10):1096.
  * Williams [2000] Williams KR.  Biomechanics in Sport; Wiley Online Library.  p. 161.
  * Williams and Cavanagh [1987] Williams KR, Cavanagh PR.  Relationship between distance running mechanics, running economy, and performance.  Journal of Applied Physiology 1987;63(3):1236–45.
  * Winter [1979] Winter DA.  A new definition of mechanical work done in human movement.  Journal of Applied Physiology 1979;46(1):79–83.


11footnotetext: corresponding author

,

# Effect of medium on fundamental interactions in gravity and condensed matter

Alexander Zhuk  Valerii Shulga

###### Abstract

Recently, it was shown that the gravitational field undergoes exponential
cutoff at large cosmological scales due to the presence of background matter.
In this article, we demonstrate that there is a close mathematical analogy
between this effect and the behavior of the magnetic field induced by a
solenoid placed in a superconductor.

Keywords:cosmology; scalar perturbations; gravitational potential; magnetic
field; superconductor.

##  1 Introduction

It seems quite natural that the presence of the medium influences the
propagation of fundamental interactions. The simplest example is the Debye
screening of the electric field of an individual particle in a plasma by
particles of opposite sign. Here, the potential produced by an external point
charge has the form of the Yukawa potential (but not the Coulomb one) with the
Debye screening length (see, e.g., [1]). A similar screening mechanism of the
electron charge due to vacuum polarization takes place in quantum
electrodynamics (see, e.g., [2]). The Anderson-Higgs mechanism is another
example of the influence of the medium on fundamental interactions, which are
carried by gauge fields. In this case, after symmetry breaking, the Higgs
vacuum field acts as a medium [3, 4, 5]. As a result of interaction with this
medium, the initially massless gauge fields gain mass [6]. It is also known
that medium in the form of the superconductor affects the electromagnetic
interaction. For example, external magnetic field undergoes the exponential
cutoff inside the superconductor due to the Meissner effect (see, e.g., [7]).

The examples above did not concern the gravitational interaction between
massive bodies. It is known that in a vacuum in the weak field limit the
gravitational potential satisfies the Poisson equation and has the form of
Newton’s potential [8]. From a naive point of view, since all masses have the
same sign and are attracted to each other, one should hardly expect a
screening of the gravitational interaction, as, for example, for electric
charges in a plasma. However, it was demonstrated recently [9, 10, 11] that
medium in the case of gravity also plays important role. It was shown that,
due to the interaction of the gravitational potential with background matter,
there is an exponential cutoff of the gravitational interaction at large
cosmological scales. In section 2 we reproduce this result. For many, this
result turned out to be rather unexpected. Therefore, in this paper, in
section 3, we present a close mathematical analogue of this phenomenon by the
example of the magnetic field induced by a solenoid placed in a
superconductor.

##  2 Screening of the gravitational interaction in cosmology

We consider the Universe containing the cosmological constant $\Lambda$ and
filled with discrete point-like gravitating sources (galaxies and the group of
galaxies) with comoving mass density

$\rho=\sum_{n}\rho_{n}=\sum_{n}m_{n}\,\delta(\mathbf{r}-\mathbf{r}_{n})\,,$
(2.1)

where $\mathbf{r}=(x^{1},x^{2},x^{3})$ is comoving distance. This is our
medium. Such matter has a dust-like equation of state and the average energy
density $\bar{\varepsilon}=\bar{\rho}c^{2}/a^{3}$ where comoving averaged mass
density $\bar{\rho}=\mathrm{const}$, $c$ is the speed of light and $a$ is the
conformal factor. The corresponding background metric is described by
Friedmann-Lema$\mathrm{\hat{\i}}$tre-Robertson-Walker (FLRW) one.

The discrete inhomogeneities perturb the FLRW metric [12, 13]:

$ds^{2}=a^{2}\big{[}(1+2\Phi)d\eta^{2}-(1-2\Phi)\delta_{\alpha\beta}\,dx^{\alpha}dx^{\beta}\,\big{]}\,,$
(2.2)

where we restrict ourselves to scalar perturbations in conformal Newtonian
gauge. Scalar function $\Phi(\eta,\mathbf{r})$ is the gravitational potential
created at the point with the radius-vector $\mathbf{r}$ by all gravitating
masses in the Universe [8]. The perturbed Einstein equations are [12, 13]:

$\Delta\Phi-3{\mathcal{H}}\left(\Phi^{\prime}+{\mathcal{H}}\Phi\right)=\frac{1}{2}\kappa
a^{2}\delta\varepsilon,$ (2.3)
$\Phi^{\prime}+{\mathcal{H}}\Phi=-\frac{1}{2}\kappa
a^{2}\bar{\varepsilon}v\,,$ (2.4)
$\Phi^{\prime\prime}+3{\mathcal{H}}\Phi^{\prime}+\left(2{\mathcal{H}}^{\prime}-{\mathcal{H}}^{2}\right)\Phi=0\,,$
(2.5)

where $\Delta\equiv\delta^{\alpha\beta}\partial_{\alpha}\partial_{\beta}$ is
the Laplace operator, the prime denotes the conformal time $\eta$ derivative,
${\mathcal{H}}\equiv(da/d\eta)/a=(a/c)H$ and $H\equiv(da/dt)/a$ is the Hubble
parameter, $v(\eta,\textbf{r})$ is the peculiar velocity potential and
$\kappa\equiv 8\pi G_{\\!N}/c^{4}$, where $G_{\\!N}$ is the gravitational
constant. The energy density fluctuation reads [14, 15]:

$\delta\varepsilon=\frac{c^{2}}{a^{3}}\delta\rho+\frac{3\bar{\rho}c^{2}}{a^{3}}\Phi\,,$
(2.6)

where $\delta\rho(\eta,\mathbf{r})\equiv\rho-\bar{\rho}$ is the fluctuation of
the mass density (2.1) around its constant average value $\bar{\rho}$.

Eq. (2.4) demonstrates that the peculiar velocities affect the gravitational
potential. If we neglect this influence (i.e.
$\Phi^{\prime}+{\mathcal{H}}\phi=0$), then equation (2.3) takes the form

$\Delta\Phi-\frac{a^{2}}{\lambda^{2}}\Phi=\frac{\kappa
c^{2}}{2a}\delta\rho\,,$ (2.7)

where the screening length

$\lambda\equiv\sqrt{\frac{2a^{3}}{3\kappa\bar{\rho}c^{2}}}\,.$ (2.8)

With the help of the transformation (to remove the $\bar{\rho}$ contribution
on the RHS of Eq. (2.7))

$\phi=c^{2}a\Phi-\frac{4\pi
G_{N}\bar{\rho}}{a^{2}}\lambda^{2}=c^{2}a\Phi-\frac{c^{2}a}{3}$ (2.9)

Eq. (2.7) is reduced to

$\Delta\phi-\frac{a^{2}}{\lambda^{2}}\phi=4\pi G_{N}\rho\,.$ (2.10)

For the mass density (2.1), we can easily solve this Helmholtz equation, and
applying transformation (2.9) obtain:

$\Phi=\frac{1}{3}-\frac{\kappa c^{2}}{8\pi
a}\sum_{n}\frac{m_{n}}{|\mathbf{r}-\mathbf{r}_{n}|}\exp\left(-\frac{a|\mathbf{r}-\mathbf{r}_{n}|}{\lambda}\right)\,.$
(2.11)

It is worth noting that the physical distance is $R=ar$. The term 1/3 (which
is due to $\bar{\rho}$ in $\delta\rho$) plays an important role since only
with this term the averaged over all volume value of the gravitational
potential $\bar{\Phi}$ is equal to zero as it should be for fluctuations [9].

In Eq. (2.11), we neglect the peculiar velocities of the inhomogeneities.
However, they also play an important role [16, 17] and must be taken into
account. For the considered model, as was shown in [16], it is sufficient in
(2.7), (2.9)-(2.11) to replace $\lambda$ with $\lambda_{\mathtt{eff}}$ and
additionally in (2.11): $1/3\to 1/3(\lambda_{\mathtt{eff}}/\lambda)^{2}$ where

$\lambda_{\mathtt{eff}}=\sqrt{\frac{c^{2}a^{2}H}{3}\int\frac{da}{a^{3}H^{3}}}\,.$
(2.12)

To get this result, we should take into consideration Eq. (2.5). This
screening length (as well as $\lambda$) depends on time. For example, for the
standard $\Lambda$CDM model at present time $\lambda_{\mathtt{eff}}=2.57$ Gpc
[16].

Therefore, the gravitational potential $\Phi$ satisfies the Helmholtz
equation, not the Poisson equation. This is due to the interaction of the
gravitational potential with the medium. We can see it directly from Eq. (2.6)
where the term $\sim\bar{\rho}\Phi$ describe this interaction. Due to the
peculiar velocity, Eq. (2.3) also acquires an additional term proportional to
$\Phi$ [16]. If the medium is absent that corresponds to the limit
$\bar{\rho}\to 0,\,v\to 0$, then the screening lengths $\lambda$ and
$\lambda_{\mathtt{eff}}$ tends to infinity, and the Yukawa potentials in
(2.11) are reduced to the Newton’s ones without screening of the gravitational
interaction.

##  3 Solenoid in a superconductor. Screening of the induced magnetic field.

In this section, in order to present the mathematical analog of the screening
effect described above, we render some of equations of the paper [18] in a
form suitable for our purpose. Following this paper, we consider a thin
solenoid placed in a superconductor. Thin means that the diameter of the
solenoid is much smaller than the magnetic field penetration length
$\lambda_{\mathtt{m}}$. It is well known that the magnetic field of the
solenoid $\mathbf{B}_{\mathtt{sol}}$ is absent from the outside it, but the
vector potential $\mathbf{A}_{\mathtt{sol}}$ is not equal to zero. The
interaction of this potential with the superconducting medium induces a
current $\mathbf{J}_{\mathtt{sc}}$, which, in turn, leads to the appearance of
an induced magnetic field $\mathbf{B}_{\mathtt{sc}}$. Thus, the Maxwell
equation has the form222In this section, we use the system of units adopted in
book [7].

$\mathtt{curl}\mathbf{B}_{\mathtt{tot}}=\mathtt{curl}\left(\mathbf{B}_{\mathtt{sc}}+\mathbf{B}_{\mathtt{sol}}\right)=\mathbf{J}_{\mathtt{sc}}+\mathbf{J}_{\mathtt{sol}}\,.$
(3.1)

Since outside the solenoid
$\mathbf{B}_{\mathtt{sol}},\mathbf{J}_{\mathtt{sol}}=0$, we get

$\mathtt{curl}\mathbf{B}_{\mathtt{sc}}=\mathbf{J}_{\mathtt{sc}}\,,$ (3.2)

where in the London limit the superconducting current density is [7, 18]

$\mathbf{J}_{\mathtt{sc}}=-\frac{1}{\lambda^{2}_{\mathtt{m}}}\left(\frac{1}{q}\nabla\theta+\mathbf{A}_{\mathtt{tot}}\right)\,.$
(3.3)

Here,
$\mathbf{A}_{\mathtt{tot}}=\mathbf{A}_{\mathtt{sc}}+\mathbf{A}_{\mathtt{sol}}$,
$\theta$ is the phase of the order parameter and the magnetic field
penetration length

$\lambda_{\mathtt{m}}=\frac{1}{q\sqrt{n_{\mathtt{s}}}}\,,$ (3.4)

whrere $n_{\mathtt{s}}$ is the superfluid density, parameter $q$ defines the
superconducting flux quanta (see, e.g., Eq. (3.7) below) and in the real
superconductor $q=2e/(\hbar c)$ [7]. The absence of a superconducting medium
corresponds to the limit $n_{\mathtt{s}}\to
0\Rightarrow\lambda_{\mathtt{m}}\to\infty$. Expression (3.4) is an analogue of
cosmological formula (2.8) (and, accordingly, formula (2.12)). In Eq. (3.3)
the term $\lambda^{-2}_{\mathtt{m}}\mathbf{A}_{\mathtt{sol}}\sim
n_{\mathtt{s}}\mathbf{A}_{\mathtt{sol}}$ describes the interaction of the
solenoid magnetic field with the superconducting medium just as the term
$\sim\bar{\rho}\Phi$ on the RHS of Eq. (2.6) describes the interactions of the
gravitational potential with the cosmological medium.

Now, applying curl operation to both sides of (3.3), we obtain

$\mathbf{B}_{\mathtt{sc}}-\lambda^{2}_{\mathtt{m}}\Delta\mathbf{B}_{\mathtt{sc}}=0\,,$
(3.5)

where we took into account that outside of the solenoid
$\mathtt{curl}\nabla\theta=0$ and $\mathbf{B}_{\mathtt{sol}}=0$. $\Delta$ is
the Laplace operator in flat space. To solve this equation, we need to define
the boundary conditions. Let the solenoid be extended along the $z$-axis.
Obviously, due to the cylindrical symmetry the induced magnetic field inside
the superconductor is also parallel to the $z$-axis:
$\mathbf{B}_{\mathtt{sc}}(r)={B}_{\mathtt{sc}}(r)\hat{z}$, where $\hat{z}$ is
the unit vector along $z$-axis. In cylindric coordinates, $\mathbf{r}$ is the
radius-vector in the $xy$-plane (it is worth noting that in the previous
section $\mathbf{r}$ denotes the comoving three-dimensional radius-vector). At
distances $r>>\lambda_{\mathtt{m}}$, the superconducting current goes to zero:
$\mathbf{J}_{\mathtt{sc}}\to 0$. Therefore, at this distances Eq. (3.3) reads

$\mathbf{A}_{\mathtt{tot}}=-\frac{1}{q}\nabla\theta\,.$ (3.6)

Integrating both sides of this equation over an area inside the contour
$r=\mbox{\rm const}\,$, and performing the Stokes area-to-contour
transformation for the RHS, we find

$\Phi_{\mathtt{tot}}=-\frac{2\pi}{q}N\equiv-\Phi_{0}N\,,\quad
N=0,1,2,\ldots\,,$ (3.7)

where $\Phi_{\mathtt{tot}}=\Phi_{\mathtt{sc}}+\Phi_{\mathtt{sol}}$ is the
total magnetic flux consisting of the sum of the magnetic fluxes of the
induced magnetic field and the magnetic field inside the solenoid. $\Phi_{0}$
is the superconducting flux quanta. Therefore,

$\Phi_{\mathtt{sc}}=\Phi_{\mathtt{tot}}-\Phi_{\mathtt{sol}}\,.$ (3.8)

This is our boundary condition. We can include it directly into Eq. (3.5):

$B_{\mathtt{sc}}-\lambda^{2}_{\mathtt{m}}\Delta
B_{\mathtt{sc}}=\Phi_{\mathtt{sc}}\delta(\mathbf{r})\,,$ (3.9)

where we took into account 2D cylindrical symmetry of the problem and,
consequently, $\Delta$ is a radial Laplace operator. Obviously, integrating
this equation over an area inside the contour $r=\mbox{\rm const}\,$ we arrive
at identity. Eq. (3.9) is the Helmholtz one (similar to Eq. (2.10)), and has
the decreasing solution

$B_{\mathtt{sc}}=\frac{\Phi_{\mathtt{sc}}}{2\pi\lambda^{2}_{\mathtt{m}}}K_{0}\left(\frac{r}{\lambda_{\mathtt{m}}}\right)\,,$
(3.10)

where $K_{0}$ is the modified Bessel function. The induced magnetic field
behaves asymptotically as follows:

$B_{\mathtt{sc}}(r\to 0)\sim-\ln(r)\,,\quad
B_{\mathtt{sc}}(r\to\infty)\sim\frac{1}{\sqrt{r}}\exp(-r/\lambda_{\mathtt{m}})\,$
(3.11)

This behavior reflects the cylindrical symmetry of the model. For example,
Yukawa’s potential has been transformed:
$(1/r)\exp(-r/\lambda_{\mathtt{m}})\to(1/\sqrt{r})\exp(-r/\lambda_{\mathtt{m}})$.
As expected, the screening length coincides with the magnetic field
penetration length $\lambda_{\mathtt{m}}$. Formula (3.11) is 2D analog of Eq.
(2.11).

##  4 Conclusion

In this paper, we have touched upon the problem of the influence of the medium
on fundamental interactions. First, on the basis of articles [9, 10, 11], we
showed that as a result of the interaction of the gravitational field with the
cosmological medium, the gravitational potential is subject to exponential
screening on large cosmological scales. Then, following the model considered
in paper [18], we have traced a close analogy between the interaction of the
gravitational field with the cosmological medium and the interaction of the
magnetic field of a solenoid with a superconducting medium. As a result of
this interaction, the induced magnetic field in the superconductor undergoes
exponential screening at distances exceeding the magnetic field penetration
length.

## Acknowledgements

The authors are grateful to Boris Svistunov for fruitful discussions and
valuable comments.

## References

  * [1] Brydges, D. C., and Martin, Ph. A. (1999). Coulomb systems at low density. J. Stat. Phys. 96, 1163–1330. arXiv:cond-mat/9904122 [cond-mat.stat-mech]. doi: 10.1023/A:1004600603161.
  * [2] Merches, I., Tatomir, D., and Lupu, R. E. (2019). Basics of Quantum Electrodynamics. NY, USA: Taylor & Francis.
  * [3] Anderson, P.W. (1963). Plasmons, gauge invariance, and mass. Phys. Rev. 130 (1), 439-442. doi:10.1103/physrev.130.439.
  * [4] Englert, F., and Brout, R. (1964). Broken symmetry and the mass of gauge vector mesons. Phys. Rev. Lett.. 13, 321-323. doi: 10.1103/PhysRevLett.13.321.
  * [5] Higgs, P. W. (1964). Broken symmetries and the masses of gauge bosons. Phys. Rev. Lett. 13, 508-509. doi: 10.1103/PhysRevLett.13.508.
  * [6] Linde, A. D. (1990). Particle Physics and Inflationary Cosmology (CONTEMPORARY CONCEPTS IN PHYSICS SERIES). NY, USA: Taylor & Francis.
  * [7] Svistunov, B., Babaev, E., and Prokof’ev N. (2015). Superfluid States of Matter. NY, USA: Taylor & Francis.
  * [8] Landau, L.D., and Lifshitz, E.M. (2000). The Classical Theory of Fields (Course of Theoretical Physics Series,V.2). Oxford, UK: Pergamon Press.
  * [9] Eingorn, M. (2016). First-order cosmological perturbations engendered by pointlike masses. Astrophys. J. 825, 84. arXiv:1509.03835 [gr-qc]. doi: 10.3847/0004-637X/825/2/84
  * [10] Eingorn, M., Kiefer, C., and Zhuk, A. (2016). Scalar and vector perturbations in a universe with discrete and continuous matter sources. J. Cosmol. Astropart. Phys. 2016, 032. arXiv:1607.03394 [gr-qc]. doi: 10.1088/1475-7516/2016/09/032
  * [11] Eingorn, M., Kiefer, C., and Zhuk, A. (2017). Cosmic screening of the gravitational interaction. Int. J. Mod. Phys. D 26, 1743012. arXiv:1711.01759 [gr-qc]. doi: 10.1142/S021827181743012X.
  * [12] Mukhanov, V.F., Feldman, H.A., and Brandenberger, R.H. (1992). Theory of cosmological perturbations. Phys. Rept. 215, 203-333. https://doi.org/10.1016/0370-1573(92)90044-Z.
  * [13] Gorbunov, D.S., and Rubakov, V.A. (2011) Introduction to the theory of the early universe: Cosmological perturbations and inflationary theory. Singapore: World Scientific.
  * [14] Eingorn, M., and Zhuk, A. (2012). Hubble flows and gravitational potentials in observable Universe. J. Cosmol. Astropart. Phys. 2012, 026. arXiv:1205.2384 [astro-ph.CO]. DOI: 10.1088/1475-7516/2012/09/026
  * [15] Eingorn, M., and Zhuk, A. (2014). Remarks on mechanical approach to observable Universe. J. Cosmol. Astropart. Phys. 2014, 024. arXiv:1309.4924 [astro-ph.CO]. doi: 10.1088/1475-7516/2014/05/024.
  * [16] Canay, E., and Eingorn, M. (2020). Duel of cosmological screening lengths. Phys. Dark Univ. 29, 100565. arXiv:2002.00437 [gr-qc]. doi: 10.1016/j.dark.2020.100565.
  * [17] Canay, E., Eingorn, M., McLaughlin II, A., Arapoğlu, A.S., and Zhuk, A. (2022). Effect of peculiar velocities of inhomogeneities on the shape of gravitational potential in spatially curved universe. arXiv:2201.07561 [gr-qc].
  * [18] Mangel, I., Kapon, I., Blau, N., Golubkov, K., Gavish, N., and Keren, A. (2020). Stiffnessometer: A magnetic-field-free superconducting stiffness meter and its application. Phys. Rev. B  102, 024502. arXiv:1705.00624 [cond-mat.supr-con]. doi: 10.1103/PhysRevB.102.024502.